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ABSTRACT 


Wireless  sensors  can  be  worn  on  soldiers  or  installed  on  vehicles  to  form 
distributed  sensor  networks  to  locate  the  source  of  sniper  fire.  A  two-step  source 
localization  process  is  proposed  for  this  sniper  detection  task.  The  time  difference  of 
arrival  (TDOA)  for  the  acoustic  signals  received  by  the  sensors  is  first  estimated  using 
the  generalized  cross  correlation  (GCC)  method.  The  estimated  TDOA  values  are  then 
used  by  the  hybrid  spherical  interpolation/maximum  likelihood  (SI/ML)  estimation 
method  to  estimate  the  shooter  location.  A  simulation  model  has  been  developed  in 
MATLAB  to  study  the  performance  of  the  hybrid  SI/ML  estimation  method.  A  wireless 
sensor  network  is  simulated  in  NS-2  to  study  the  network  throughput,  delay  and  jitter. 
Simulation  results  indicate  that  the  estimation  accuracy  can  be  increased  by  increasing 
the  number  of  sensors  or  the  inter-sensor  spacing.  The  constraint  of  small  inter-sensor 
spacing  on  wearable  sensors  is  found  to  degrade  the  estimation  accuracy,  but  vehicular 
configuration  providing  larger  inter-sensor  spacing  can  help  improve  the  estimation 
accuracy.  The  sensor  topology  should  be  well  represented  in  all  three  dimensions  to 
obtain  desired  estimation  accuracy.  The  estimation  accuracy  is  not  adversely  affected  by 
sensor  node  failures  or  location  perturbations.  The  NS-2  simulation  results  indicate  that 
the  wireless  sensor  network  has  low  delay  and  can  support  fast  information  exchange 
needed  in  counter-sniper  applications. 
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EXECUTIVE  SUMMARY 


In  urban  warfare,  sensor  networks  can  be  effectively  used  to  detect,  track,  and 
monitor  targets  of  interest.  Sniper  localization  through  the  collaborative  effort  of  multiple 
wireless  acoustic  sensors  is  of  interest  in  this  study.  Wireless  sensors  worn  on  soldiers  or 
installed  on  vehicles  can  be  deployed  to  form  a  distributed  sensor  array  to  locate  the 
source  of  sniper  fire.  A  two-step  source  localization  process  is  proposed  in  this  work  for 
counter-sniper  applications:  the  time  difference  of  arrival  (TDOA)  values  are  first 
determined  using  the  generalized  cross  correlation  (GCC)  method  and  then  the  TDOA 
values  are  used  by  a  hybrid  spherical  interpolation/maximum  likelihood  (SI/ML) 
estimation  method  to  determine  the  shooter  location.  A  simulation  model  has  been 
developed  in  MATLAB  to  study  the  performance  of  the  hybrid  SI/ML  estimation 
method.  A  wireless  sensor  network  is  simulated  in  NS-2  to  study  the  network  throughput, 
delay  and  jitter  of  IEEE  WPAN  and  WEAN  networks. 

The  effects  of  number  of  sensor  nodes,  inter-sensor  spacing,  topology,  noise, 
source-sensor  distance,  sensor  node  failure,  location  perturbation  on  the  location 
estimation  performance  have  been  studied  through  MATEAB  simulation.  Simulation 
results  show  that  the  error  variance  can  be  reduced  by  increasing  the  number  of  sensor 
nodes  or  the  inter-sensor  spacing.  To  increase  the  location  estimation  accuracy,  the 
number  of  sensors  as  well  as  the  inter-sensor  spacing  should  be  increased. 

Eor  the  case  of  wearable  sensors,  the  constraint  of  small  inter-sensor  spacing  on 
the  body  degrades  the  estimation  performance.  On  the  other  hand,  vehicular 
configuration  providing  larger  vertical  and  horizontal  inter-sensor  spacing  can  help 
improve  the  location  estimation  accuracy.  The  sensor  node  topology  should  be  well 
represented  in  all  three  dimensions  to  obtain  desired  location  estimation  performance. 

The  source  location  estimation  accuracy  is  not  adversely  affected  by  sensor  node 
failures  or  location  perturbations.  Results  show  that  sensor  node  failures  lead  to  a  gradual 
degradation  of  accuracy  and  not  a  sudden  total  system  failure.  This  demonstrates  the 
advantage  of  using  distributed  wireless  sensor  networks  to  provide  greater  system 


xvii 


robustness.  It  is  also  observed  that  location  perturbations  do  not  adversely  affect  the 
location  estimation.  This  illustrates  the  flexibility  of  the  network  configuration  and 
demonstrates  the  suitability  for  potential  field  deployment. 

Network  simulation  results  based  on  the  NS-2  package  show  that  the  wireless 
sensor  network  can  support  a  throughput  between  20  kbps  and  100  kbps,  a  low  average 
delay  of  less  than  10  ms  and  a  low  average  jitter  in  the  range  of  5  ms  to  10  ms.  These  are 
suitable  for  fast  response  network-based  sniper  detection  applications  that  require  low 
delay  and  fast  data  exchange. 
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logical  link  control 

LR-WPAN 

Low  Rate  Wireless  Personal  Area  Network 

LS 

least-squares 

MAC 

medium  access  layer 

ML 

maximum  likelihood 

OSI 

open  systems  interconnection 

PHAT 

phase  transform 

PHY 

physical  layer 

RD 

range  difference 

RDOA 

range  difference  of  arrival 

RMS 

root  mean  square 

RSS 

received  signal  strength 

SI 

spherical  interpolation 

SRP 

Steered  Response  Power 

SRP-PHAT 

Steered  Response  Power-Phase  Transform 

SSCS 

service  specific  convergence  sub-layer 

TOA 

time  of  arrival 

TDOA 

time  difference  of  arrival 

WLAN 

wireless  local  area  network 

WPAN 

wireless  personal  area  network 

error  mean 

error  variance 

r 

RMS  error 
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I. 


INTRODUCTION 


Low-cost,  low-power,  multi-functional  wireless  sensor  networks  are  becoming 
readily  available.  The  sensor  nodes  are  small  in  size  and  are  able  to  sense,  process  data, 
and  communicate  with  each  other  wirelessly. 

Wireless  sensor  networks  have  several  potential  applications,  especially  in  the 
area  of  military  operations  [1][2].  In  urban  warfare,  sensors  can  be  deployed  in  buildings 
that  have  been  cleared  to  prevent  reoccupation.  Sensors  can  be  used  to  detect  physical 
phenomena,  such  as  motion,  acceleration,  acoustic  and  light,  and  the  collected  and 
processed  sensor  data  can  be  transmitted  through  wireless  links  and  displayed  on  mobile 
devices  carried  by  soldiers. 

In  urban  warfare,  sensor  networks  can  be  effectively  used  to  detect,  track,  and 
monitor  targets  of  interest.  Sniper  localization  through  the  collaborative  effort  of  multiple 
wireless  acoustic  sensors  is  of  interest  in  this  study.  An  array  of  acoustic  sensors, 
microphones  or  piezoelectric  pressure  sensors,  are  deployed  to  detect  and  process 
supersonic  shockwaves  and  muzzle  blast  produced  by  the  gun  shot  [3].  The  measured 
acoustic  data  are  further  processed  using  source  localization  techniques  in  order  to 
estimate  the  location  of  the  shooter. 

A.  THESIS  OBJECTIVE 

This  thesis  focused  on  source  localization  using  wireless  sensor  networks. 
Wireless  sensors  worn  on  soldiers  or  installed  on  vehicles  can  be  deployed  to  form  a 
distributed  sensor  array  to  locate  the  source  of  sniper  fire.  To  improve  the  accuracy  of  the 
source  localization,  a  large  number  of  sensors  can  be  deployed  to  cover  a  large  area  of 
interest  [4].  As  illustrated  in  Figure  1,  a  distributed  wireless  sensor  network,  forming  a 
large  array  of  acoustic  sensors,  is  used  to  detect  and  locate  sniper  fire. 

A  two-step  source  localization  process  is  proposed  in  this  work  for  counter-sniper 
applications:  the  time  difference  of  arrival  (TDOA)  values  are  first  determined  using  the 
generalized  cross  correlation  (GCC)  method  and  then  the  TDOA  values  are  used  by  a 
hybrid  spherical  interpolation/maximum  likelihood  (SI/ML)  estimation  method  to 
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determine  the  shooter  location.  A  simulation  model  is  developed  in  MATLAB  to  study 
the  performance  of  the  hybrid  SI/ML  estimation  method.  A  wireless  sensor  network  is 
simulated  in  NS-2  to  study  the  network  throughput,  delay  and  jitter  of  IEEE  WPAN  and 
WEAN  networks. 
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Eigure  1 .  Sniper  detection  and  localization  using  a  distributed  wireless  sensor 

network. 


B.  RELATED  WORK 

Several  physical  phenomena  can  be  used  for  sniper  detection  purposes.  Pauli  et  al. 
have  built  the  VIPER  system  that  utilizes  an  infrared  camera  to  detect  the  muzzle  flash  of 
the  weapon  [5].  It  is  augmented  with  a  microphone  to  detect  the  muzzle  blast  for  range 
estimation.  Both  sensors  require  direct  line  of  sight. 

Commercial  acoustic  sniper  detection  systems,  such  as  BBN’s  Bullet  Ears  and 
Metravib’s  PIEAR,  are  also  based  on  measurement  of  the  time  of  arrival  (TOA)  of 
muzzle  blasts  and  shockwaves  produced  by  the  gun  shot.  BBN’s  Bullet  Ears  system 
utilizes  one  or  two  small  arrays  of  microphones,  providing  estimates  of  range  to  the 
shooter  and  the  caliber,  speed  and  trajectory  of  the  projectile  [6]. 

Simon  et  al.  have  developed  an  experimental  network-based  counter-sniper 
system  called  PinPtr  [4].  The  system  utilizes  an  ad  hoc  wireless  sensor  network  built 
from  inexpensive  sensor  nodes.  After  deployment,  the  sensor  nodes  synchronize  their 
clocks,  perform  self-localization  and  wait  for  acoustic  events.  The  sensors  can  detect 
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muzzle  blasts  and  acoustic  shockwaves  and  measure  their  time  of  arrival  (TOA). 
Utilizing  a  message  routing  service,  the  TOA  measurements  are  delivered  to  a  laptop 
computer,  where  the  sensor  signal  processing  algorithm  calculates  the  shooter  location 
estimate.  The  proposed  solution  utilizes  time  of  arrival  data  of  the  measured  shockwaves 
and  muzzle  blasts.  From  the  measurements  and  the  sensor  positions  a  four-dimensional 
consistency  function  is  formed.  A  quick  search  algorithm  finds  the  maximum  of  this 
function.  The  location  corresponding  to  the  maximum  is  the  shooter  position  estimate. 


C.  THESIS  OUTLINE 

The  organization  of  the  thesis  is  as  follows:  chapter  II  introduces  the  wireless 
sensor  network  as  well  as  the  various  source  localization  techniques  and  sniper  detection 
approaches.  Chapter  III  describes  the  concept  of  source  localization  using  wireless 
sensor  networks.  The  two-step  source  localization  technique  proposed  for  counter-sniper 
applications  will  be  discussed.  Chapter  IV  presents  the  simulation  model  of  the  two-step 
source  localization  process,  followed  by  the  simulation  results  of  the  hybrid  SI/ML 
estimation  method  using  the  MATLAB  package  simulating  a  sniper  detection  scenario. 
The  throughput,  delay  and  jitter  of  the  wireless  sensor  network  are  simulated  using  the 
NS-2  network  simulator.  Chapter  V  concludes  the  thesis  and  highlights  future  work  for 
further  investigation.  Appendix  includes  the  MATLAB  source  codes  used  in  the 
simulation  studies. 
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II.  INTRODUCTION  TO  WIRELESS  SENSOR  NETWORK, 
SOURCE  LOCALIZATION  AND  SNIPER  DETECTION 


This  chapter  discusses  the  various  aspects  of  distributed  wireless  sensor  networks, 
acoustic  source  localization  and  current  sniper  detection  systems.  We  will  introduce  the 
concept  of  using  wireless  sensor  networks  to  perform  source  localization  for  the  purpose 
of  sniper  detection. 

A.  DISTRIBUTED  WIRELESS  SENSOR  NETWORKS 

Low-cost,  low-power,  multi-functional  wireless  sensor  networks  are  becoming 
readily  available.  The  wireless  sensor  nodes  are  small  in  size  and  are  able  to 
communicate  with  each  other  wirelessly  by  forming  an  ad  hoc  network. 

Wireless  sensor  networks  are  self-organizing  networks  with  short-range 
communication  and  multi- hop  routing.  Their  topologies  frequently  change  due  to 
mobility,  fading  and  node  failures.  There  are  usually  limitations  on  transmit  power, 
memory,  and  computational  capability  in  a  sensor  node  [1]. 

The  performance  of  a  wireless  sensor  network  generally  depends  on  the  following 
factors  [1]:  energy  efficiency,  latency,  accuracy,  fault  tolerance,  scalability, 
synchronization  and  localization,  and  throughput.  Due  to  the  interdependence  of  energy 
consumption,  delay,  and  throughput,  all  these  performance  issues  and  metrics  are  tightly 
coupled. 

For  this  thesis  research,  the  emphasis  is  not  on  any  of  the  above  performance 
metrics.  Instead,  we  will  focus  on  the  issue  of  distributed  array  processing  for 
determining  the  location  of  a  signal  source  of  interest. 

I.  Architecture 

Wireless  sensor  networks  can  be  organized  in  two  architectures:  flat  or  tiered.  In  a 
flat  architecture,  all  sensor  nodes  are  peers  and  are  homogeneous  in  form  and  function. 
However,  in  a  tiered  architecture  [7],  as  shown  in  Figure  2,  sensor  nodes  form  a  hierarchy 
in  which  a  sensor  node  at  a  given  level  performs  a  specific  set  of  tasks  on  behalf  of  a 
subset  of  sensor  nodes  in  the  level  below. 
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Tier 


Figure  2.  Tiered  Architecture  in  Low  Rate  Wireless  Personal  Area  Networks  (LR- 

WPAN)  (from  [7]). 

In  practice,  very  few  sensor  networks  are  entirely  flat.  Typically,  a  sensor  network 
will  connect  to  an  internet  via  gateway  nodes.  In  addition,  sensor  networks  are  often  not 
physically  homogeneous.  For  example,  a  network  may  become  heterogeneous  due  to 
uneven  battery  drain  across  nodes.  Tiered  architectures  are  usually  employed  to  take 
advantage  of  unevenly  distributed  resources  by  assigning  resource-intensive  roles  to 
resource-rich  nodes. 

2.  Standards 

Efforts  are  under  way  to  standardize  the  various  layers  of  wireless  sensor  network 
communication  protocols.  IEEE  802.15.4  is  the  standard  defined  for  Eow  Rate  Wireless 
Personal  Area  Networks  (ER-WPAN)  [8],  with  the  ZigBee  Alliance  [9]  as  its  marketing 
and  compliance  certification  organization. 

A  ER-WPAN  is  a  simple,  low-cost  communication  network  that  allows  wireless 
connectivity  in  applications  with  limited  power  and  relaxed  throughput  requirements.  The 
main  objectives  of  a  ER-WPAN  are  ease  of  installation,  reliable  data  transfer,  short-range 
operation,  low  cost,  and  a  reasonable  battery  life,  while  maintaining  a  simple  and  flexible 
protocol  suite. 
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a.  LR-WPAN  Architecture 

The  LR-WPAN  architecture  is  specified  in  terms  of  a  number  of  layers 
based  on  the  open  systems  interconnection  (OSI)  seven-layer  model,  as  shown  in  Figure 
3.  Each  layer  is  responsible  for  one  part  of  the  standard  and  offers  services  to  higher 
layers. 


Figure  3.  IEEE  802.15.4  ER-WPAN  device  architecture  (from  [8]). 

The  IEEE  802.15.4  standard  defines  the  physical  layer  (PHY)  and  the 
medium  access  layer  (MAC).  An  IEEE  802.2  Type  I  logical  link  control  (EEC)  can 
access  the  MAC  layer  through  the  service  specific  convergence  sub-layer  (SSCS).  The 
PHY  layer  specifies  an  operating  frequency  of  2.4  GHz  with  a  basic  bit  rate  of  250  kbps. 
There  are  alternate  PHY  specifications  for  915  MHz  and  868  MHz  that  operate  at  lower 
data  rates  of  40  kbps  and  20  kbps,  respectively.  Details  of  the  ER-WPAN  modulation  and 
spreading  formats  are  summarized  in  Table  1. 
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PHY(MHz) 

Erequency 

band 

(MHz) 

Spreading  parameters 

Data  parameters 

Chip  rate 
(kchips/s) 

Modulation 

Bit  rate 
(kb/s) 

Symbol 

rate 

(ksym/s) 

Symbols 

868 

868-868.6 

300 

BPSK 

20 

20 

Binary 

915 

902-928 

600 

BPSK 

40 

40 

Binary 

2450 

2400- 

2483.5 

2000 

0-QPSK 

250 

62.5 

16-ary 

Orthogonal 

Table  1. 


IEEE  802.15.4  ER-WPAN  frequency  bands  and  data  rates  (from  [8]). 


b.  LR-  WPAN  Network  T apologies 

Depending  on  the  application  requirements,  the  ER-WPAN  may  operate  in 
one  of  the  two  topologies  shown  in  Eigure  4:  Star  Topology  or  Peer-to-Peer  Topology.  In 
the  Star  Topology,  the  communication  is  established  between  devices  and  a  single  central 
controller,  called  the  Personal  Area  Network  (PAN)  coordinator.  Peer-to-Peer  Topology 
also  has  a  PAN  coordinator;  however,  it  differs  from  the  star  topology  in  that  any  device 
can  communicate  with  any  other  device  as  long  as  they  are  within  range  of  one  another. 

Peer-to-Peer  Topology  allows  more  complex  network  formations  to  be 
implemented,  such  as  mesh  networking.  A  peer-to-peer  network  can  be  ad  hoc,  self¬ 
organizing  and  self-healing.  It  may  also  allow  multiple  hops  to  route  messages  from  any 
device  to  any  other  device  on  the  network. 
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Figure  4.  Star  and  Peer-to-Peer  topologies  for  a  IEEE  802. 15.4  ER-WPAN  (from 

[8]). 


To  form  a  larger  network,  supporting  hundreds  and  thousands  of  sensor 
nodes,  a  single  sensor  cluster  can  reorganize  and  form  a  mesh  of  multiple  neighboring 
clusters.  Once  predetermined  application  or  network  requirements  are  met,  the  PAN 
coordinator  may  instruct  a  device  to  become  the  cluster  head  of  a  new  cluster  adjacent  to 
the  first  one.  Other  devices  gradually  connect  and  form  a  multi-cluster  network  structure, 
such  as  the  one  seen  in  Eigure  5.  The  advantage  of  a  multi-cluster  structure  is  increased 
coverage  area;  the  disadvantage  is  an  increase  in  message  latency. 


Eigure  5.  IEEE  802.15.4  ER-WPAN  cluster  tree  network  (from  [8]). 
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B.  ACOUSTIC  SOURCE  UOCAUIZATION 


Sound  localization  can  be  accomplished  by  using  observed  differences  in  the 
sound  signals  reeeived  at  different  observation  points,  sueh  as  angle  of  arrival  (AOA), 
time  of  arrival  (TOA),  time  difference  of  arrival  (TDOA),  received  signal  strength  (RSS) 
or  steered  response  power  (SRP),  to  estimate  the  direetion  and  actual  location  of  the 
sound  source.  Several  methods,  usually  based  on  SRP  and  TDOA,  have  been  developed 
for  estimating  an  acoustical  source  location. 

Algorithms,  like  Steered  Response  Power- Phase  Transform  (SRP-PHAT)  [10], 
are  robust  in  the  presence  of  reverberation.  However,  SRP-PHAT  can  be  quite  complex 
requiring  the  calculation  of  a  large  number  of  test  points  in  the  region  of  possible  source 
locations.  The  location  is  chosen  to  be  the  point  that  produces  the  highest  steered 
response  power. 

Alternatively,  the  problem  can  be  implemented  as  a  two-step  localization  process 
[11].  In  a  two-step  localization  approach,  the  time  difference  of  arrival  (TDOA)  is  first 
determined  by  performing  cross-correlation  on  the  received  signals  [12].  The  TDOA 
values  are  then  used  to  estimate  the  loeation. 

There  are  generally  three  classes  of  TDOA  estimators  [12],  including  the  general 
cross-correlation  (GCC)  approaeh,  the  maximum  likelihood  (ML)  approaeh,  and  the 
phase  transform  (PHAT)  approach.  All  these  approaches  attempt  to  calculate  the  cross¬ 
correlation  in  an  optimal  or  suboptimal  manner,  and  then  seleet  the  time  index  of  the  peak 
of  the  result  to  be  the  TDOA  estimate. 

The  two-step  localization  algorithms  are  quite  fast,  but  they  lack  robustness  [12]. 
Frequently,  errors  do  occur  in  the  time  delay  estimates  due  to  reflections  of  the  sound 
source,  whieh  are  sometimes  greater  in  energy  than  the  direct  signal.  The  direct  path  can 
be  obstrueted  or  attenuated  beeause  of  source  and  microphone  direetivity.  Erroneous  time 
delay  estimates  eause  large  errors  in  location  estimation. 

For  the  counter-sniper  application  studied  in  this  thesis,  we  propose  the  use  of  the 
two-step  localization  process  due  to  its  speed  and  the  limited  computational  capabilities 
of  the  sensors. 
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C.  SNIPER  DETECTION 

A  gun  shot  produces  several  detectable  phenomena,  sueh  as  muzzle  blast  and 
shoekwaves  [3]  [6]. 

Muzzle  blast  is  the  aeoustic  signature  associated  with  the  ejection  of  the  bullet 
from  the  sniper’s  rifle.  The  muzzle  blast  is  a  loud,  eharaeteristie  noise  originating  from 
the  end  of  the  muzzle  and  propagating  spherieally  away  at  the  speed  of  sound.  The 
muzzle  blast  can  be  detected  with  acoustie  sensors  at  ranges  from  several  hundred  meters 
out  to  more  than  a  kilometer. 

Typical  rifles  fire  projectiles  at  supersonic  velocities,  thereby  produeing  acoustie 
shoeks  along  their  trajeetory.  Shoekwaves  ean  be  deteeted  aeoustically  at  ranges  from 
hundreds  of  meters  out  to  more  than  a  kilometer  and  can  be  used  to  accurately  determine 
projectile  trajectories.  The  shock  waveform  is  distinctive  and  cannot  be  produced  by  any 
other  natural  phenomenon. 

As  illustrated  in  Figure  6,  the  muzzle  blast  and  shockwaves  can  be  deteeted  and 
processed  by  distributed  wireless  sensor  networks  equipped  with  acoustie  sensors. 
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Figure  6.  Shock  wave  front  and  muzzle  wave  detection  by  a  group  of  distributed 

wireless  sensors  (after  [4]). 
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As  illustrated  in  Figure  7,  wireless  sensors  may  be  mounted  on  vehieles  for  easy 
and  fast  deployment.  The  wireless  sensors  ean  also  be  worn  on  soldiers  or  integrated  into 
their  helmets. 


Figure  7.  Wireless  network-based  eounter-sniper  system  consisting  of  vehicular  and 

wearable  sensors. 

We  have  introduced  the  basics  of  wireless  sensor  networks,  described  the  concept 
of  source  localization  for  acoustic  signal  sources  and  proposed  the  concept  of  using 
wireless  sensor  networks  to  detect  sniper  fire.  In  the  next  chapter,  we  will  discuss  in 
greater  detail  the  various  source  localization  techniques.  We  will  discuss  the  two-step 
source  localization  technique  proposed  for  counter-sniper  applications. 
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III.  SOURCE  LOCALIZATION  USING  WIRELESS  SENSOR 

NETWORKS 


This  chapter  presents  the  details  of  souree  loealization  using  distributed  sensor 
networks.  We  will  discuss  the  two-step  souree  localization  technique  proposed  for 
counter-sniper  applications. 

Source  loealization  is  the  process  of  determining  the  spatial  loeation  of  a  source 
based  on  multiple  observations  of  the  emitted  signal  source.  The  primary  goal  of  source 
loealization  is  accuracy.  The  aeeuraey  of  loeation  estimates  depends  on  a  number  of 
faetors,  ineluding  quantity  and  quality  of  sensors  used,  sensor  placement  relative  to  eaeh 
other  and  the  signal  sources  to  be  analyzed,  ambient  noise  and  reverberation,  and  number 
of  aetive  sources  and  their  speetral  eontent.  The  performance  of  loealization  techniques 
generally  improves  with  the  number  of  sensors  in  the  array,  especially  when  adverse 
environmental  eonditions  are  present. 

A.  TWO-STEP  LOCALIZATION  PROCESS 

Figure  8  shows  the  schematic  diagram  of  the  two-step  source  localization 
technique  considered  in  this  thesis.  The  generalized  cross  correlation  (GCC)  teehnique  is 
used  to  estimate  the  time  differenee  of  arrival  a  pair  of  received  signals  m^(t)  and 

m2(0 . 

The  fj2  values  are  then  eonverted  into  range  difference  of  arrival  d^2  values  using 
the  relationship 

d^2  —  '^12^ 

where  i>  =  345m/s  is  the  speed  of  sound  [12].  Figure  9  shows  the  range  differenee  of 
arrival  at  two  aeoustie  sensors.  The  Jj2  values  along  with  the  knowledge  of  the  sensor 

positions  are  then  used  to  generate  hyperbolic  curves  whieh  are  then  intersected  to  obtain 
a  souree  loeation  estimate. 
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Figure  8.  Two-step  source  localization  model  using  GCC  and  hybrid  SI/ML 
methods.  Step  1  performs  the  TDOA  estimation  using  GCC;  Step  2  performs  the 
location  estimation  in  two  stages:  estimate  determined  by  SI  method  in  stage  1  is 
used  as  initial  values  for  ML  method  in  stage  2. 


Acoustic  waves 


Figure  9.  Range  difference  of  arrival  at  two  acoustic  sensors. 


Location  estimation  can  be  performed  using  either  a  closed-form  least- 
squares  (LS)  [13]  or  a  maximum  likelihood  (ML)  method  [14].  The  LS  method  is  fast, 
and  a  closed-form  solution  can  be  obtained  in  a  single  iteration.  The  spherical 
interpolation  (SI)  method  [13]  is  a  widely  used  LS  implementation  in  location  estimation 
problems.  In  contrast,  the  ML  method  requires  iterative  gradient-descent  search,  and  a 
good  initial  guess  close  to  the  solution  is  required  in  order  to  avoid  divergence. 
However,  the  ML  estimator  is  asymptotically  unbiased  and  provides  a  more  accurate 
estimate  at  low  noise  levels. 
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In  this  work,  we  propose  the  use  of  a  hybrid  SI/ML  algorithm  [14]  for 
location  estimation  using  the  two-stage  process,  as  illustrated  in  Figure  8:  First,  we 
calculate  the  source  location  using  spherical  interpolation  (SI)  method.  The  estimate 
of  the  source  location  is  then  used  as  initial  value  for  the  calculation  of  source 
location  5^^  using  the  maximum  likelihood  (ML)  method.  A  hybrid  SI/ML  algorithm  has 

been  proposed  by  Ajdler  et  al.  [14]  to  perform  acoustic  source  localization  in  a  small 
room  of  size  3mx5mx  3  m.  In  this  thesis,  a  similar  hybrid  SFML  algorithm  is  applied 
for  acoustic  source  localization  in  an  outdoor  environment  instead  to  detect  sniper  fire. 

At  low  noise  levels,  the  SI  method  provides  a  rough  estimate  quickly.  The 
ML  method  is  computationally  demanding  and  requires  the  use  of  an  iterative  technique, 
such  as  Gauss-Newton  and  Levenberg-Marquardt  [15].  The  convergence  speed  of  the  ML 
method  can  be  increased  by  using  the  rough  estimate  from  the  SI  method  as  initial 
estimate.  Additionally,  in  situations  where  the  ML  method  fails  to  converge,  especially  at 
high  noise  levels,  the  SI  method  can  still  provide  a  reasonable  estimate. 

Other  closed-form  estimation  methods  like  spherical  intersection  (SX) 
[13],  planar  intersection  (PX)  [13],  and  clustering-based  techniques  [16]  can  also  be  used 
to  produce  the  initial  guess  of  the  source  location.  However,  the  SI  method  has  greater 
noise  immunity  than  the  SX  method  and  the  SI  method  has  lower  error  variance  but 
slightly  higher  bias  than  the  PX  method.  Clustering-based  techniques  are  computationally 
less  demanding  than  the  SI  method  but  are  only  applicable  to  sources  that  are  very  far 
away. 

Alternatively,  the  ML  optimization  problem  can  be  solved  by  using 
alternating  projection  method  [17]  to  avoid  exhaustive  multi-dimensional  search. 
However,  this  is  computationally  more  intensive  than  the  closed-form  methods  and  there 
is  still  no  guarantee  that  the  convergence  to  a  global  minimum  can  be  obtained. 
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1.  Time  Difference  of  Arrival  Estimation 

Let  m,  (t)  and  mj  (t)  represent  the  time-delayed  signals  of  a  source  signal  received 
at  two  sensors  as  shown  in  Figure  9.  If  Mj(ty)  and  (<»)  represent  the  Fourier 
transforms  of  received  signals  m^(t)  and  m^it) ,  respectively,  the  TDOA  Tcan  be 
estimated  using  the  GCC  given  as  [12] 

f  =  argmaxf  W {co)M^[co)Ml[co)e^'^^ dO) 

P  J-oo 

where  T  is  an  estimate  of  the  delay  between  original  source  signal  and  the  two  sensors, 
and  VF(<»)is  the  GCC  weighting  function. 

Common  GCC  weighting  functions  include  maximum  likelihood  (ML)  and  phase 
transform  (PHAT)  [12]: 

I  (<»)p  (<»)p  -I- 1^2  (<»)p  |Mi  (<»)p 

^PHAT  {^)  ~  I  /  \  *  /  O 

Mj  {(o)M^  {o))\ 

where  Ni((o)  and  AiCm)  are  the  estimated  noise  spectral  density  for  the  first  and  second 
sensors,  respectively. 

Phase  transform  (PHAT)  uses  only  the  phase  information  derived  from  their 
respective  signal  power  spectral  density  [12].  This  phase-only  procedure  has  the 
whitening  effect  on  the  signals  and  emphasizes  primarily  the  channel  effect.  The  PHAT 
has  a  peak  at  the  relative  time  delay  of  the  two  signals  and  is  used  to  estimate  the  TDOA 
information. 

The  ML  weights  require  knowledge  about  the  spectrum  of  the  sensor-dependent 
noises.  The  PHAT  weights  do  not  require  this  knowledge,  and  hence  has  been  employed 
more  often  due  to  its  simplicity.  We  choose  the  GCC-PHAT  method  due  to  its  ability  to 
reduce  the  degradation  due  to  reverberation  [10]. 
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2.  Least-squares  Spherical  Interpolation  (SI)  Method 

The  Spherical  Interpolation  method  [13]  is  chosen  for  the  closed- form  least- 
squares  (LS)  source  location  estimation  based  on  the  measured  range  difference  of  arrival 
(RDOA)  values.  Figure  10  shows  the  notations  used  in  SI  method  for  a  single  source 
and  two  sensors  i  and  the  reference  sensor  1 . 

Let  N  denote  the  number  of  sensors.  The  vector  of  the  spatial  coordinates  for  the 
i-th  sensor  is  denoted  by^.  =  [x.  y.  zj^and  the  estimated  position  of  the  source  is 

denoted  by^  =  [v^  Sensor  1  is  arbitrary  selected  as  the  reference  node  and 

thus  the  position  of  sensor  1  is  set  as  the  origin  of  the  coordinate  system  (i.e.,  ^  =  0) 
The  distance  between  the  source  and  sensor  i  is  denoted  by  the  Euclidean  norm 


and  the  distance  from  the  origin  to  the  point  v,  is  denoted  by  R.  =  .  Similarly,  the 


estimated  distance  from  the  origin  to  the  source  is  denoted  by  =  |  ^  | .  The  measured 
RDOA  between  sensor  i  and  sensor  1  is  denoted  by 


fori  =  2,  ,  N . 
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Figure  10.  Notations  used  in  Spherieal  Interpolation  (SI)  method  for  a  single  source, 

sensor  i  and  the  reference  sensor  1  (from  [13]). 


Assuming  that  sensor  1  is  selected  as  the  reference  node  (i.e.,  ^  =  0 ),  it  can  be 
shown  [13]  that  the  source  location  estimate  is  given  as 


1-2 


'j(s^ P,^WP,‘-sY'  PfWP,^S. 


1  1  •  •  •  T  dd^ 

where  the  projection  matnx  =  1  - 

d  d 


S  = 


pd.  _ 

1^2  U21 

RI  -  dl^ 


pd.  _  Jd. 

^(Vl. 


the  range  difference  of  arrival  vector,  d  = 


d 
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d. 
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^2  3^2  ^2 


the  sensor  location  matrix,  S 

y3 

^3 

.  9 

yN 

_ 

••• 

and  the  weight  matrix,  W  = 

>^21 

^22 

^2N 

_Wm 

Wn2 

^NN 

The  range  estimate  is  given  as  the  Euclidean  norm  from  the  origin  at  sensor  1  to 
the  source 

K  =  ||1||  = 

and  the  bearing  estimate  cosQ^  is  defined  as  the  direction  cosines  from  the  origin  at 
sensor  1  to  the  source 

A 

cosi2  =^  =  -^. 

K  HI 

The  bearing  estimate  cos  is  equivalent  to  the  unit  vector  of  the  source  location 
estimate  s_ ,  as  shown  in  Figure  1 1 . 


source 


Figure  1 1 .  Range  and  bearing  estimate  of  source  location. 
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3.  Maximum  Likelihood  Estimator 

The  errors  r]  in  the  measurement  of  the  range  difference  of  arrival  (RDOA)  d-j 

are  assumed  to  be  zero  mean  Gaussian  and  independent  for  each  sensor.  We  define  E  as 
the  covariance  matrix  of  the  noise  vector,  A  as  the  vector  of  all  the  exactly  measured 
RDOAs,  andr  as  the  vector  containing  the  noisy  measurements: 

r  =  A+7. 

Assuming  that  sensor  1  is  chosen  as  the  reference  node,  it  can  be  shown  [14]  that 
the  likelihood  function  of  F given  the  source  location^ is 


/(r,i) 


1 


where  the  RDOA  error  vector  h  (^)  is  represented  by 


(1.1) 


the  noise  covariance  matrix  is  represented  by 


and  the  ML  cost  function  is  the  exponent 


(1.2) 


(1.3) 


(1.4) 


From  Equation  (1.1)  and  Equation  (1.4),  it  can  be  seen  that  maximizing  the 
likelihood  function /(F,^)  is  equivalent  to  minimizing  the  ME  cost  function  7^^  (^) . 

Hence,  the  source  location  estimate  is  obtained  by  minimizing  the  ME  cost  function 
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4L=argmin  . 

£ 

The  cost  function  7^^  (^ )  can  be  minimized  using  standard  numerical 

optimization  methods.  The  Levenberg-Marquardt  method  [15]  provided  by  MATLAB  is 
used  to  solve  the  minimization  problem. 


4.  Cramer  Rao  Lower  Bound 

The  minimum  mean- squared  error  for  any  estimate  of  a  non-random  parameter  is 
given  by  the  Cramer  Rao  Bound  [18].  The  Cramer  Rao  Bound  (CRB)  gives  a  lower 
bound  on  the  error  variance  of  any  unbiased  estimate.  The  variance  of  any  unbiased 
estimator  of  ^is  bounded  by 


where  F{£)  is  the  Fischer  information  matrix,  expressed  as 
F{i)  =  E 


(1.5) 


piog/(r,5) 

~diogf(r,sy 

T  ' 

(1.6) 


where  Pis  the  vector  containing  the  noisy  RDOA  measurements  and  /(F,^)  is  the 
likelihood  function  of  F  given  the  source  location  5  . 

From  Equation  (1.1)  and  Equation  (1.6),  it  can  be  shown  that  [14] 


FU) 


dh[s) 

T 

y-1 

dh[^) 

(1.7) 


where  J  (^y)  is  the  Jacobian  matrix 


J(i) 
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From  Equation  (1.5)  and  Equation  (1.7),  it  can  be  shown  that  the  minimum 
varianee  of  the  estimates  of  the  souree  position  is  equal  to  the  traee  of  the  inverse  of  the 
Eiseher  information  matrix  [14] 


minE 


(l-iXI-if 


=  tr(F-'U)). 


(1.8) 


Eaeh  diagonal  element  of  the  inverse  Eiseher  information  matrix  corresponds  to 
the  minimal  varianee  of  one  eoordinate  of  the  souree  position.  The  total  varianee  on  the 
estimation  of  the  souree  position  is  the  trace  of  the  inverse  Eiseher  information  matrix.  It 
ean  be  shown  [14]  that  the  bound  on  the  varianee  of  the  souree  loeation  estimate  is 
direetly  proportional  to  the  noise  varianee. 

This  ehapter  diseussed  the  two-step  souree  loealization  teehnique  proposed  for 
our  eounter-sniper  applieations.  The  generalized  eross  eorrelation  seheme  used  to 
determine  the  time  differenee  of  arrival  was  briefly  diseussed  and  the  least-squares 
spherieal  interpolation  (SI)  and  maximum  likelihood  (ME)  teehniques  have  been 
summarized.  In  the  next  chapter,  we  will  present  the  simulation  results  of  the  hybrid 
SI/ME  estimation  method  as  well  as  study  the  network  throughput,  delay  and  jitters  of  the 
wireless  networks. 
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IV.  SIMULATION  RESULTS  AND  ANALYSIS 


In  this  study,  the  two-step  source  localization  process,  as  discussed  in  Chapter  III, 
is  used  to  estimate  the  sniper  location  through  the  use  of  wireless  sensor  nodes  forming  a 
distributed  sensor  network  to  detect  sniper  fire.  The  generalized  cross  correlation  (GCC) 
technique,  as  shown  in  Figure  8,  is  first  used  to  estimate  the  time  difference  of  arrival 
(TDOA)  of  the  received  acoustic  signals  at  the  sensors.  The  TDOA  values,  as  shown  in 
Figure  9,  are  then  converted  into  range  difference  of  arrival  (RDOA)  values  to  be  used  in 
the  hybrid  SI/ML  estimation  method  to  estimate  the  source  location. 

This  chapter  presents  the  simulation  model  of  the  two-step  source  localization 
process,  followed  by  the  simulation  results  of  the  hybrid  SFML  estimation  method  using 
the  MATLAB  package  simulating  a  sniper  detection  scenario.  The  MATLAB  code  used 
for  the  simulation  is  provided  in  Appendix.  The  throughput,  delay  and  jitter  of  the 
wireless  sensor  network  are  simulated  using  the  NS-2  network  simulator. 

A.  SIMULATION  MODEL 

A  simulation  model  is  developed  in  MATLAB  to  study  the  two-step  localization 
process.  Figure  12  shows  the  simulation  block  diagram  for  the  study  of  the  two-step 
localization  process.  The  study  of  the  TDOA  estimation  using  the  GCC  technique  is 
carried  out  by  examining  the  received  signals  at  the  sensors.  The  study  of  the  hybrid 
SI/ML  estimation  is  carried  out  separately  by  representing  TDOA  values  as  RDOA  and 
modeling  RDOA  estimation  errors  as  Gaussian  noise  [14].  The  acoustic  source  is  not 
explicitly  simulated  in  this  work. 


Figure  12.  Simulation  block  diagram  for  the  study  of  the  two-step  localization 

process  (based  on  Figure  8  in  Chapter  III). 
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1.  Time  Difference  of  Arrival  Estimation 

In  the  simulation  model,  based  on  the  known  positions  of  the  source  and  sensors, 
time-delayed  copies  of  a  source  signal  are  generated  for  each  sensor.  We  now  consider 
the  scenario  of  two  sensors  (see  Figure  9).  Figure  13  shows  the  simulated  received 
signals  at  sensors  1  and  2.  The  received  signal  at  sensor  2  is  the  same  as  the  signal 
received  at  sensor  1  but  is  time  delayed  by  =  40  ms. 

The  TDOA  between  the  two  sensors  is  determined  using  the  GCC  method  as  has 
been  discussed  in  Chapter  III.  By  applying  the  GCC  technique  on  the  two  received 
signals,  the  TDOA  estimate  Tjj  is  determined  to  be  40.6  ms,  which  is  close  to  the  actual 

value  of  Tjj  =40  ms. 


Figure  13.  Received  signals  at  two  sensors.  The  received  signal  at  sensor  2  is  same  as 
the  signal  received  at  sensor  1  but  is  time  delayed  by  =  40  ms. 

2.  Hybrid  SI/ML  Location  Estimation 

To  study  the  hybrid  SFML  location  estimation  method  proposed  in  the  two-step 
localization  process  for  counter-sniper  applications,  the  TDOA  information  T-. 

determined  using  the  GCC  technique  is  now  converted  into  RDOA  values  d..  using  the 

relationship  d^-  where  V  =  345  m/s  is  the  speed  of  sound  [12].  In  the  simulation 
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study,  however,  the  GCC  method  is  not  used.  Instead,  we  have  chosen  to  generate  the 
measured  RDOA  values  using  an  additive  white  Gaussian  model  as  described  below. 

The  true  RDOA  dij  values  can  be  determined  using  the  known  source  and  sensor 
positions  that  are  defined  in  the  simulation  scenario.  Figure  14  shows  the  range  difference 
of  arrival  for  a  single  source  and  two  sensors;  sensor  1  is  chosen  to  be  the  reference 
sensor.  The  locations  of  the  source,  sensor  2,  and  sensor  1  are  (vj,  >'2>^2) 

(0,0,0),  respectively.  The  distances  Djand  are  the  Euclidean  norms  between  the 
source  and  sensor  1  and  the  source  and  sensor  2,  respectively.  The  RDOA  (ij2  between 
sensor  1  and  sensor  2  is  given  by 

d^2  —  —  D2  —  '^[2^ 

where  t^2  *^he  magnitude  of  the  difference  between  the  time  of  arrival  at 

sensor  1  and  the  time  of  arrival  ^2  at  sensor  2. 


^12  |h  h 


Figure  14.  Range  difference  of  arrival  d^2  for  a  single  source,  sensor  2  and  reference 
sensor  1  located  at  (^2,  y2’^2)  (0,0,0) ,  respectively. 


The  measured  RDOA  values  d^j  are  modeled  as 

d.j=d.j  +  n 

where  n  is  a  white  Gaussian  noise  with  zero  mean  and  variance  al .  Figure  15  illustrates 
the  additive  white  Gaussian  noise  model  used  to  generate  the  measured  RDOA  values 
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djj .  The  true  RDOA  values  d^j  are  determined  using  the  known  sensor  and  source 
positions  as  illustrated  in  Figure  14. 


0 


n 


Figure  15.  Modeling  of  range  difference  of  arrival  (RDOA). 


The  d.j  values  are  then  used  in  the  hybrid  SI/ML  estimation  method  to  estimate 


the  source  location.  As  discussed  in  Chapter  III,  the  hybrid  SFML  estimation  method 
consists  of  two  stages.  As  shown  in  Figure  12,  the  first  stage  of  the  hybrid  SFML  method 
produces  the  location  estimate  using  the  SI  method.  The  SI  estimate  is  then  used  to 
provide  initial  values  for  the  iterative  ML  method  used  in  the  second  stage.  In  the 
simulation  study,  location  estimation  performance  of  both  SI  and  ML  are  studied  using 
the  performance  metrics  described  below. 

Monte  Carlo  simulation  runs  are  carried  out  and  the  estimated  source  location 
coordinates  s-  are  compared  with  the  actual  source  location  coordinates  s  to  obtain  the 
performance  metrics  of  error  mean 


error  variance 
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and  root  mean  square  (RMS)  error  given  by  square  root  of  the  second  moment 


y=^jU^  +  C7^ 

where  M  is  the  number  of  simulation  runs. 

Figure  16  shows  location  estimates  of  100  Monte  Carlo  simulation  runs  of  the  SI 
and  ML  methods  in  a  scenario  with  the  source  at  100  m  from  the  origin;  the  sensors  are 
spaced  10  m  apart.  Table  2  shows  the  error  mean  // ,  error  variance  (7^  and  RMS  error  y 
for  24  and  72  sensor  nodes  for  the  SI  and  ML  methods.  The  notations  Xs_SI,  Ys_SI, 
Zs_SI,  Rs_SI  and  Bearing_SI  refer  to  the  individual  X,  Y,  Z  range  and  bearing  readings, 
respectively,  for  the  SI  method  used  in  the  first  stage  of  the  hybrid  SI/ML  method. 
Similarly,  the  notations  Xs_ML,  Ys_ML,  Zs_ML,  Rs_ML  and  Bearing_ML  refer  to  the 
individual  X,  Y,  Z  range  and  bearing  readings,  respectively,  for  the  ML  method  used  in 
the  second  stage.  From  Figure  16,  the  estimated  locations  do  not  coincide  with  the  actual 
source  location  at  (0,100,0) .  For  the  case  of  24  sensor  nodes,  the  error  variance  is 
large  and  the  position  estimates  are  widely  scattered.  For  the  case  of  72  sensor  nodes, 
however,  the  ML  method  produces  estimation  with  a  low  RMS  error  of  less  than  3  m,  and 
the  position  estimates  are  closer  to  the  source  at  (0,100,0) . 
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Figure  16.  Simulation  output  of  the  location  estimation  using  the  SI  and  ML  methods 
in  a  seenario  with  the  source  at  100  m  from  the  origin  and  the  sensors  spaced  10 
m  apart:  (a)  24  sensor  nodes  and  (b)  72  sensor  nodes.  The  number  of  Monte  Carlo 

runs  is  100. 
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Performance  Metrics 

Number  of  Nodes 

24 

72 

Error  Mean 

in  meters 

Xs SI 

-1.021 

-7.6879 

Ys Sl 

-2.0383 

20.688 

Zs Sl 

-1.3407 

-0.7423 

Rs SI 

-1.8955 

20.992 

Bearing SI 

2.8112 

3.7754 

Error  mean 

in  meters 

Xs ML 

0.007947 

-0.76399 

Ys ML 

-0.37713 

1.6088 

Zs ML 

-0.01764 

-0.17491 

Rs ML 

-0.32973 

1.6182 

Bearing ML 

1.535 

0.54189 

Error  variance 
in  square 

meters 

Xs SI 

11.914 

11.224 

Ys SI 

359.54 

10.753 

Zs SI 

7.5425 

4.1336 

Rs SI 

353.59 

12.112 

Bearing SI 

11.599 

1.9837 

Error  variance 
in  square 

meters 

Xs ML 

5.4121 

1.2701 

Ys ML 

196.71 

5.7233 

Zs ML 

2.828 

0.11484 

Rs ML 

195.4 

5.7962 

Bearing ML 

1.4399 

0.29564 

RMS  error 

Y^i  in  meters 

Xs SI 

3.5996 

8.3862 

Ys Sl 

19.071 

20.946 

Zs SI 

3.0561 

2.1644 

Rs SI 

18.899 

21.279 

Bearing SI 

4.4161 

4.0296 

RMS  error 
in  meters 

Xs ML 

2.3264 

1.3615 

Ys ML 

14.03 

2.883 

Zs ML 

1.6818 

0.38136 

Rs ML 

13.983 

2.9008 

Bearing_ML 

1.9484 

0.76765 

Table  2.  Mean  ji ,  variance  ,  and  RMS  value  J  of  tho  position  errors  for  the  SI  and  ML 
methods  in  a  scenario  with  the  source  at  100  m  from  the  origin  and  the  sensors 
spaced  10  m  apart.  The  number  of  Monte  Carlo  runs  is  100. 


29 


B.  HYBRID  SI/ML  ESTIMATION  METHOD  SIMULATION 


The  effects  of  the  number  of  sensor  nodes,  inter-sensor  spacing,  topology,  noise, 
source-sensor  distance,  sensor  node  failure,  and  location  perturbation  on  the  estimation 
performance  of  the  hybrid  implementation  of  the  closed-form  least-squares  spherical 
interpolation  (SI)  method  and  the  iterative  maximum  likelihood  (ML)  method  are  studied. 
The  main  objective  of  this  portion  of  the  simulation  is  to  study  suitable  wearable  and 
vehicular  sensor  node  configurations  for  the  counter-sniper  application  as  well  as  to 
determine  the  required  number  of  sensors  and  the  inter-sensor  spacing  to  accurately 
estimate  a  source  of  sniper  fire  located  a  few  hundred  meters  away. 

Figure  17  shows  the  schematic  of  a  sensor  node  topology  using  the  three- 
dimensional  XYZ  Cartesian  coordinate  system.  By  setting  different  source  and  sensor 
positions,  different  simulation  scenarios  can  be  constructed.  Note  that  for  the  sensor 
node  topology  with  regular  inter-sensor  spacing,  a  simplified  notation  of  KxLxM  is 
used  to  describe  a  topology  that  has  K ,  L  and  M  sensor  nodes  along  the 
X,Y  and  Z  plane,  respectively. 


source 


•  l(0,y,z) 


(0,0,z)  X 


(0,0,0) 


Figure  17.  Schematic  of  the  sensor  topology  using  the  three-dimensional  XYZ 


Cartesian  coordinate  system.  White  and  grey  circles  represent  sensors  that  are 
placed  in  the  two-dimensional  XY  and  the  three-dimensional  XYZ  planes, 

respectively. 
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1.  Effect  of  Number  of  Sensor  Nodes 

The  effect  of  the  number  of  sensor  nodes  on  the  location  estimation  performance 
is  studied  using  the  scenario  of  a  source  located  100  m  from  the  origin;  the  sensors  are 
spaced  10  m  apart.  Several  topologies  are  formed  using  9  (i.e.,  3x1x3),  12  (i.e.,  2x2x3), 
24  (i.e.,  4x2x3),  36  (i.e.,  6x2x3),  48  (i.e.,  8x2x3),  60  (i.e.,  10x2x3)  and  72  (i.e., 
12x2x3)  sensor  nodes.  A  3x1x3  topology  represents  a  single  row  of  three  sensor 
columns  with  three  sensors  each;  a  12x2x3  topology  represents  two  rows  of  twelve 
sensor  columns  with  three  sensors  each. 

Figure  18  and  Table  3  show  the  error  mean  //  ,  error  variance  (J^and  RMS  error 
^for  different  sensor  node  topologies.  For  the  cases  of  less  than  12  sensor  nodes,  the 
errors  of  the  ML  were  very  large  and  are  not  shown  in  the  plots  in  Figure  18. 

It  can  be  seen  that  as  the  number  of  sensor  nodes  is  increased,  the  error  variance 
decreases  for  both  SI  and  ML  methods  used  in  the  two  stages  of  the  hybrid  SI/ML 
algorithm.  For  example,  when  the  number  of  sensors  increases  from  24  to  48,  the  error 
variance  a  for  the  Rs  values  in  the  SI  stage  decreases  quickly  from  354  m  to  36  m 
(approximately).  The  error  variance  <7^  of  the  ML  stage  is  lower  than  the  first  SI  stage 
and  decreases  from  195  m  to  7  m  (approximately).  It  is  observed  that  when  the  error 
variance  of  the  SI  method  in  the  first  stage  is  high  (i.e.,  for  the  case  of  9  and  12  sensor 
nodes),  it  led  to  highly  inaccurate  results  for  the  ML  method  in  the  second  stage. 

The  effect  of  the  bias  introduced  by  the  SI  method  is  clearly  seen  when  the 
number  of  sensors  is  very  high  (i.e.,  48  sensors  or  more).  However,  the  ML  method  in  the 
second  stage  is  tolerant  to  the  bias  introduced  by  the  SI  stage.  For  example,  the  error 
mean //for  Rs  values  using  the  SI  method  increased  from  9  m  to  21  m  (approximately) 
when  the  number  of  sensors  increased  from  48  to  72;  however,  the  RMS  error  y  for  Rs 
values  for  the  ML  method  remained  below  3.1  m. 
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Figure  18.  Effect  of  number  of  sensor  nodes  on  the  location  estimation  performance 
of  the  hybrid  SEME  method:  (a)  error  mean  fi ,  (b)  error  variance  ,  and  (c) 
RMS  error  j .  The  readings  for  9  and  12  sensors  are  not  shown  due  to  highly 
erroneous  results.  Results  are  based  on  100  Monte  Carlo  runs. 
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Number 
of  Nodes 

9 

12 

24 

36 

48 

60 

72 

Error 

mean 

in  meters 

Xs Sl 

-0.99945 

-0.04349 

-1.021 

-1.9582 

-2.8292 

-4.7036 

-7.6879 

Ys SI 

-26.868 

-30.456 

-2.0383 

5.0957 

8.7197 

13.989 

20.688 

Zs Sl 

1.4638 

1.8704 

-1.3407 

-1.6174 

-1.3639 

-1.1141 

-0.7423 

Rs SI 

-23.144 

-26.308 

-1.8955 

5.2157 

8.8395 

14.155 

20.992 

Bearing 

SI 

21.487 

22.94 

2.8112 

2.3622 

2.2773 

2.7098 

3.7754 

Error 

mean 

in  meters 

Xs ML 

-3.49E-r08 

1.28E-I-09 

0.007947 

-0.77733 

-0.52786 

-0.5967 

-0.76399 

Ys ML 

1.60E-I-08 

5.49E-r08 

-0.37713 

2.9459 

1.5521 

1.4334 

1.6088 

Zs ML 

-1.91E-I-08 

7.18E-I-08 

-0.01764 

-0.34894 

-0.18401 

-0.1626 

-0.17491 

Rs ML 

4.37E-I-08 

1.57E-I-09 

-0.32973 

2.9664 

1.5599 

1.4398 

1.6182 

Bearing 

ML 

10.497 

23.621 

1.535 

0.97759 

0.54588 

0.46911 

0.54189 

Error 
variance 
ah  in 
square 

meters 

Xs SI 

10.163 

7.3159 

11.914 

12.586 

12.01 

11.316 

11.224 

Ys SI 

2088.6 

2156.9 

359.54 

75.292 

34.693 

18.157 

10.753 

Zs Sl 

22.946 

23.819 

7.5425 

6.8466 

5.3038 

4.4692 

4.1336 

Rs Sl 

1560.7 

1591.6 

353.59 

75.748 

35.639 

19.304 

12.112 

Bearing 

si 

1970.4 

2098.1 

11.599 

1.8786 

1.8306 

1.9483 

1.9837 

Error 

variance 

in 

square 

meters 

Xs ML 

1.52E-I-18 

4.12E-I-16 

5.4121 

2.6975 

1.0198 

0.80249 

1.2701 

Ys ML 

2.76E-I-17 

7.71E-I-17 

196.71 

32.131 

6.962 

4.2166 

5.7233 

Zs ML 

4.79E-I-17 

7.32E-I-17 

2.828 

0.90724 

0.30799 

0.15477 

0.11484 

Rs ML 

2.26E-I-18 

1.54E-I-18 

195.4 

32.201 

6.9973 

4.2598 

5.7962 

Bearing 

ML 

447.78 

233.55 

1.4399 

0.34273 

0.20057 

0.17955 

0.29564 

RMS  error 
Ysi  in 
meters 

Xs SI 

3.341 

2.7051 

3.5996 

4.0523 

4.4737 

5.7827 

8.3862 

Ys Sl 

53.014 

55.538 

19.071 

10.063 

10.523 

14.624 

20.946 

Zs SI 

5.0089 

5.2266 

3.0561 

3.0761 

2.6766 

2.3897 

2.1644 

Rs SI 

45.786 

47.788 

18.899 

10.146 

10.667 

14.821 

21.279 

Bearing 

S1 

49.316 

51.228 

4.4161 

2.7311 

2.6489 

3.0482 

4.0296 

RMS  error 
Yml  in 
meters 

Xs ML 

1.28E-I-09 

2.12E-I-08 

2.3264 

1.8171 

1.1395 

1.0764 

1.3615 

Ys ML 

5.49E-r08 

9.20E-r08 

14.03 

6.3883 

3.0612 

2.5043 

2.883 

Zs ML 

7.18E-I-08 

8.98E-I-08 

1.6818 

1.0144 

0.58468 

0.42569 

0.38136 

Rs ML 

1.57E-I-09 

1.30E-r09 

13.983 

6.4032 

3.0709 

2.5165 

2.9008 

Bearing 

_ML 

23.621 

17.475 

1.9484 

1.1395 

0.70609 

0.63215 

0.76765 

Table  3.  Effect  of  number  of  sensor  nodes  on  the  location  estimation  performance  of  the 
hybrid  SI/ML  method:  error  mean  // ,  error  variance  ,  and  RMS  error  y .  The 

number  of  Monte  Carlo  runs  is  100. 
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Another  scenario  is  simulated  to  represent  a  number  of  soldiers  wearing  two 
sensors  spaced  1  m  apart.  The  sensors  on  the  soldier  are  assumed  to  be  mounted  one 
each  on  the  foot  and  the  helmet.  Each  soldier  is  stationed  10  m  apart  and  the  source  is  at 
100  m  from  the  origin.  This  leads  to  a  3-D  sensor  node  arrangement,  and  the  scenario 
considered  the  sensor  node  arrangements  in  the  range  from  4x2x2  to  12x2x2.  A 
4x2x2  indicates  a  total  eight  soldiers  wearing  two  sensors  each  in  two  rows  of  four 
soldiers  each. 

Figure  19  shows  the  error  mean,  error  variance  and  RMS  error  for  the  different 
number  of  soldiers.  It  can  be  seen  that  the  estimation  accuracy  in  the  Z  direction  is  poorer 
than  in  the  X  and  Y  directions.  For  example,  for  the  case  of  5x2x2,  the  RMS  errors  for 
the  Xs,  Ys,  and  Zs  readings  at  the  MF  stage  are  3  m,  10  m  and  14  m,  respectively.  The 
poorer  accuracy  in  the  Z  direction  is  due  to  the  constraint  that  the  wearable  sensors  are 
only  1  m  apart  in  the  Z  direction. 

It  is  observed  that  by  deploying  more  sensor  nodes,  the  accuracy  of  the  estimation 
can  be  improved.  For  example,  twelve  soldiers  with  two  sensors  each  (i.e.,  6x2x2)  can 
detect  a  sniper  100  m  away  with  an  RMS  error  of  3  m,  8  m  and  15  m  for  the  Xs,  Ys, 
and  Zs  readings,  respectively,  at  the  MF  stage.  However,  when  the  number  of  soldiers  is 
greater  than  eighteen,  the  larger  bias  introduced  by  the  SI  method  in  the  first  stage  may 
lead  to  less  accurate  estimates  for  the  MF  stage,  especially  for  the  Zs  values.  For 
example,  forty-eight  soldiers  (12x2x2)  can  detect  a  sniper  100  m  away  with  an  RMS 
error  7  of  3  m,  4  m  and  17  m  for  the  Xs,  Ys,  and  Zs  readings,  respectively,  using  the  MF 
method. 
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Figure  19.  Performance  of  the  hybrid  SI/ML  method  for  simulated  wearable 

configuration  (two  sensors  per  soldier)  with  different  number  of  soldiers:  (a)  error 
mean  // ,  (b)  error  variance  ,  and  (c)  RMS  error  y .  Results  are  based  on  100 

Monte  Carlo  runs. 


2.  Effect  of  Inter- Sensor  Spacing 

To  study  the  effect  of  inter-sensor  spacing  on  the  location  estimation 
performance,  a  scenario  consisting  of  24  sensor  nodes  (6x2x2)  is  simulated.  The  source 
is  situated  100  m  from  the  origin  and  the  inter-sensor  spacing  is  varied  from  10  m  to  50 
m.  Figure  20  and  Table  4  show  the  error  mean  ji ,  error  variance  and  RMS  error  y 
for  different  inter-sensor  spacing. 

It  is  observed  that  when  the  inter-sensor  spacing  increases,  the  estimation 
accuracy  increases.  The  error  variance  decreases  rapidly  when  the  inter-sensor 
spacing  is  increased.  For  example,  at  the  inter-sensor  spacing  of  10  m,  the  error  variance 

a  of  the  Rs  values  obtained  at  the  ML  stage  was  195  m  ;  at  the  inter-sensor  spacing  of 

2  2 

20  m,  the  error  variance  o  quickly  decreased  to  7.6  m  . 
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hybrid  SI/ML  method  :  (a)  error  mean  // ,  (b)  error  variance  <7^ ,  and  (c)  RMS 
error  y .  Results  are  based  on  100  Monte  Carlo  runs. 
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Inter¬ 

sensor 

Spacing 

1  m 

5  m 

10  m 

20  m 

30  m 

40  m 

50  m 

Mean 

error 

in  meters 

Xs SI 

1.4874 

3.8611 

-1.021 

-3.9662 

-5.0268 

-6.0904 

-7.0891 

Ys Sl 

-96.052 

-56.341 

-2.0383 

6.8829 

5.3421 

4.125 

2.4439 

Zs SI 

1.014 

2.3276 

-1.3407 

-3.667 

-5.3771 

-7.7525 

-10.657 

Rs SI 

-92.607 

-49.787 

-1.8955 

7.0553 

5.6146 

4.5986 

3.2445 

Bearing_ 

SI 

75.399 

42.603 

2.8112 

2.8987 

3.9867 

5.4073 

7.1224 

Mean 

error 

in 

meters 

Xs ML 

-6.54E-I-07 

4.01E-I-07 

0.00795 

-1.0804 

-0.6564 

-0.30166 

-0.04878 

Ys ML 

l.llE-i-08 

2.03E-r08 

-0.37713 

3.4316 

1.5148 

0.55075 

0.055739 

Zs ML 

-2.82E-I-07 

3.03E-I-06 

-0.01764 

-0.75852 

-0.50311 

-0.24626 

-0.04254 

Rs ML 

1.35E-r08 

2.46E-I-08 

-0.32973 

3.4495 

1.5218 

0.55344 

0.056847 

Bearing_ 

ML 

26.764 

10.791 

1.535 

0.91882 

0.57074 

0.34763 

0.2187 

Error 
variance 
(7gi  in 
square 

meters 

Xs SI 

0.054935 

1.19E-I-01 

11.914 

5.0661 

2.1157 

0.99307 

0.451 

Ys Sl 

244.13 

1.97E-I-03 

359.54 

14.286 

1.6202 

0.1256 

0.014497 

Zs Sl 

0.024743 

5.78E-I-00 

7.5425 

3.6117 

1.5424 

0.7841 

0.41151 

Rs SI 

208.03 

1.38E-r03 

353.59 

15.188 

1.9432 

0.22131 

0.002205 

Bearing_ 

SI 

3060.8 

3297.6 

11.599 

1.9415 

0.96781 

0.4869 

0.25618 

Error 

variance 

in 

square 

meters 

Xs ML 

1.64E-I-16 

7.94E-I-16 

5.4121 

1.2418 

0.44208 

0.21931 

0.11579 

Ys ML 

4.77E-I-16 

2.39E-I-17 

196.71 

7.6353 

0.97013 

0.20748 

0.032234 

Zs ML 

4.07E-I-15 

1.57E-I-16 

2.828 

0.81316 

0.32197 

0.17749 

0.10389 

Rs ML 

6.74E-I-16 

3.16E-I-17 

195.4 

7.7062 

0.98182 

0.2098 

0.032316 

Bearing_ 

ML 

709.58 

174.95 

1.4399 

0.27527 

0.12701 

0.055 

0.025056 

RMS 

error 

in  meters 

Xs SI 

1.64E-I-16 

5.18E-r00 

3.5996 

4.5604 

5.233 

6.1714 

7.1209 

Ys SI 

4.77E-I-16 

7.17E-I-01 

19.071 

7.8524 

5.4917 

4.1402 

2.4469 

Zs SI 

4.07E-I-15 

3.35E-I-00 

3.0561 

4.1302 

5.5186 

7.8029 

10.676 

Rs SI 

6.74E-I-16 

6.21E-I-01 

18.899 

8.0601 

5.7851 

4.6226 

3.2449 

Bearing_ 

SI 

709.58 

71.503 

4.4161 

3.2162 

4.1063 

5.4522 

7.1404 

RMS 

error 

Yml  in 
meters 

Xs ML 

1.44E-I-08 

2.85E-r08 

2.3264 

1.5521 

0.93432 

0.55705 

0.34377 

Ys ML 

2.45E-I-08 

5.29E-r08 

14.03 

4.4058 

1.8069 

0.71471 

0.18799 

Zs ML 

6.98E-I-07 

1.25E-I-08 

1.6818 

1.1783 

0.75834 

0.48799 

0.32511 

Rs ML 

2.93E-r08 

6.14E-I-08 

13.983 

4.4278 

1.816 

0.7184 

0.18854 

Bearing_ 

ML 

37.761 

17.07 

1.9484 

1.0581 

0.67287 

0.41934 

0.26998 

Table  4.  Effect  of  inter-sensor  spacing  on  location  estimation  performance  of  the  hybrid 
SEME  method:  error  mean  // ,  error  variance  ,  and  RMS  error  j .  Results  are 

based  on  100  Monte  Carlo  runs. 
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Additional  simulations  are  carried  out  to  take  into  account  the  limited  vertical 
inter-sensor  spacing  when  the  sensors  are  worn  on  soldiers  or  installed  on  vehicles.  It  is 
assumed  that  there  is  limited  height  of  less  than  2  m  for  the  soldiers  and  less  than  7  m  for 
vehicles  with  elevated  masts. 


a.  Wearable  Configuration 

To  study  the  wearable  configuration  for  the  soldiers,  another  scenario  is 
simulated  to  represent  twelve  soldiers  (6x2x2)  with  two  sensors  per  soldier  for  different 
Z-spacing  values  of  1.0  m,  1.25  m,  1.5  m,  1.75  m  and  2  m.  Figure  21  shows  the  RMS 
error  y  for  the  different  Z-spacing  values. 

It  can  be  seen  that  the  RMS  error  decreases  when  the  Z-spacing  is 
increased.  For  example,  the  RMS  errors  for  Zs  readings  at  the  ML  stage  decreased  from 
15  m  to  8  m  when  the  Z-spacing  was  increased  from  1  m  to  2  m;  however,  the  RMS 
errors  for  Xs  and  Ys  readings  decreased  slightly  from  3  m  and  8  m  to  2  m  and  7  m, 
respectively.  Hence,  it  is  recommended  that  the  two  sensors  worn  on  the  soldiers  should 
be  separated  as  far  apart  as  practically  possible;  however,  considering  the  constraint  of 
the  average  height  of  a  person  ranging  from  1.5  m  to  2  m,  the  assumption  of  1  m  for  the 
Z-spacing  separation,  say  between  the  head  to  the  waist  level,  is  a  reasonable  assumption 
for  a  wearable  configuration. 
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Figure  21.  RMS  error  7 for  different  values  of  Z-spacing  for  6x2x2  sensor  nodes, 
representing  a  wearable  configuration.  Results  are  obtained  by  averaging  output 

from  100  Monte  Carlo  simulations. 


b.  Vehicular  Configuration 

A  vehicular  configuration  can  provide  a  larger  vertical  inter-sensor 
spacing  compared  to  the  wearable  configuration.  The  second  sensor  can  be  mounted  on 
an  elevated  mast  to  provide  a  larger  sensor  spacing  in  the  Z  direction  and  will  improve 
the  estimation  accuracy. 

To  study  the  vehicular  configuration,  the  scenario  representing  eight 
vehicles  (4x2x2)  with  two  sensors  installed  vertically  on  each  vehicle  is  simulated.  The 
vehicles  are  spaced  10  m  apart  and  the  source  is  situated  100  m  from  the  origin.  The  Z- 
spacing  values  are  varied  from  1  m  to  7  m.  Figure  22  shows  the  RMS  error  for  different 
Z-spacing  values. 

When  the  Z-spacing  is  increased,  the  estimation  accuracy  increases.  For 
example,  with  a  Z-spacing  of  3  m,  the  RMS  error  y  at  the  ML  stage  was  3  m,  14  m  and  4 
m  for  Xs,  Ys,  and  Zs  readings,  respectively;  with  a  Z-spacing  of  7  m,  the  RMS  errors 
decreased  to  2  m,  12  m  and  2  m  for  X5,  Ys,  and  Zs  readings,  respectively. 
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Figure  22.  RMS  error  y  for  different  values  of  Z-spacing,  from  1  m  to  7  m, 

representing  a  vehicular  configuration.  Results  are  obtained  by  averaging  output 
from  100  Monte  Carlo  simulations. 


A  similar  scenario  consisting  of  eight  vehicles  is  repeated,  but  the 
simulation  is  carried  out  with  different  values  of  inter- vehicle  spacing,  from  10  m  to  40 
m.  Figure  23  shows  the  RMS  error  for  different  vehicular  spacing  values. 

It  can  be  seen  that  the  estimation  accuracy  can  be  improved  by  increasing 
the  inter- vehicle  spacing  beyond  10  m.  For  example,  for  an  inter- vehicle  spacing  of  10 
m,  the  RMS  error  at  the  ML  stage  was  15  m  for  the  Rs  readings;  for  an  inter- vehicle 
spacing  of  40  m,  the  RMS  error  was  significantly  reduced  to  4  m  for  the  same  reading. 

In  this  scenario,  however,  wireless  links  that  can  cover  a  longer  physical 
range  will  be  needed.  For  instance,  the  IEEE  802.11  WEAN  network  with  an  effective 
range  of  up  to  100  m  can  be  used  instead  of  the  IEEE  802.15.4  WPAN  that  has  an 
effective  range  of  up  to  only  10  m. 
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Figure  23.  RMS  error  ^'for  different  values  of  inter- vehicular  spacing  of  10  m,  20  m, 
30  m,  and  40  m,  representing  a  vehicular  configuration.  Results  are  obtained  by 
averaging  output  from  100  Monte  Carlo  simulations. 


3.  Effect  of  Topology 

After  studying  the  effect  of  the  number  of  sensor  nodes  and  the  inter-sensor 
spacing,  further  simulations  are  carried  out  to  find  out  a  suitable  sensor  topology  for  the 
counter-sniper  applications  and  the  required  number  of  sensors  to  be  worn  on  soldiers.  To 
study  the  effect  of  topology  on  the  location  estimation  performance,  different  two- 
dimensional  (4x2x1  and  8x2x1)  and  three-dimensional  topologies  (4x2x2  and 
4x2x3)  are  simulated.  Figure  24  shows  four  different  simulated  topologies  used  in  the 
study.  The  4x2x1  and  8x2x1  topologies  refer  to  the  layout  of  eight  and  sixteen  sensors 
on  the  XY  plane,  respectively.  The  4x2x2  and  4x2x3  topologies  refer  to  the  layout  of 
eight  sensor  columns  with  two  sensors  and  three  sensors,  respectively,  along  the  Z  axis. 
Each  sensor  is  separated  by  a  distance  of  10  m. 
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Figure  24.  Different  sensor  topologies:  (a)  4x2x1,  (b)  8x2x1,  (c)  4x2x2  and  (d) 
4x2x3.  White  and  grey  circles  represent  sensors  that  are  placed  in  the  two- 
dimensional  XY  and  three-dimensional  XYZ  planes,  respectively. 


Figure  25  shows  the  RMS  errors  for  the  four  sensor  topologies  considered. 
Simulation  results  indicate  that  although  the  estimation  errors  are  low  in  the  X  direction 
for  the  case  of  eight  sensors  lying  on  the  XY  plane  (i.e.,  4x2x  1  topology),  the  estimation 
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errors  in  the  Y  and  Z  directions  are  relatively  large.  For  example,  RMS  error  of  the  Xs,  Ys 
and  Zs  values  at  the  ML  stage  for  the  4x2x1  topology  was  10  m,  40  m  and  30  m, 
respectively. 
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Figure  25.  RMS  error  y  for  different  two-dimensional  (4x2x1  and  8x2x1)  and 
three-dimensional  topologies  (4x2x2  and  4x2x3).  Results  are  obtained  by 
averaging  output  from  100  Monte  Carlo  simulations. 

The  same  observation  of  large  estimation  errors  in  the  Y  and  Z  directions  can  be 
seen  when  the  number  of  sensors  lying  on  the  XY  plane  is  increased  from  eight  to  sixteen 
(i.e.,  8x2x1  topology).  For  example,  the  RMS  error  for  the  Xs,  Ys  and  Zs  values  at  the 
ML  stage  for  the  8x2x  1  topology  was  15  m,  35  m  and  30  m,  respectively.  However,  it  is 
observed  that  when  the  number  of  sensor  nodes  in  the  Z  direction  is  increased  (i.e., 
4x2x2  and  4x2x3  topologies),  there  is  significant  improvement  in  the  estimation 
accuracy  in  the  all  three  directions.  For  example,  the  4x2x1  and  8x2x1  topologies 
produce  the  largest  error  of  30  m  in  the  Zs  readings  at  the  ML  stage  while  the  4x2x2  and 
4x2x3  topologies  produce  the  smallest  error  of  3  m  in  the  Zs  readings.  Hence,  it  is 
recommended  that  the  sensor  topology  should  be  well  represented  in  all  three  directions 
to  obtain  good  estimation  accuracy  in  all  three  dimensions. 

Another  scenario  is  carried  out  to  represent  topologies  for  sensors  worn  on 
soldiers.  For  wearable  sensors,  it  is  assumed  that  there  is  a  limited  height  of  less  than  2  m 
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where  the  sensors  can  be  worn  on  the  soldier’s  foot  and  helmet.  The  simulation  is  carried 
out  by  having  one  topology  of  two  wearable  sensors  spaced  1  m  to  1.5  m  apart  and 
another  topology  of  three  wearable  sensors  spaced  0.5  m  to  0.75  m  apart. 


Figure  26  shows  the  RMS  error  ^'for  different  number  of  vertical  sensors  per 


soldier:  two  wearable  sensors  spaced  1  m  to  1.5  m  apart  (6x2x2)  and  three  wearable 
sensors  spaced  0.5  m  to  0.75  m  apart  (6x2x3).  It  can  be  seen  that  the  estimation 
accuracy  is  comparable  for  the  two-sensor  (6x2x2)  and  three-sensor  wearable 
configurations  (6x2x3).  For  example,  the  topologies  of  6x2x2  (Z-spacing  =  1  m)  and 
6x2x3  (Z-spacing  =  0.5  m)  produced  similar  RMS  errors  of  3  m,  8  m  and  15  m  for  X5, 
Ys  and  Zs  readings,  respectively,  at  the  ML  stage.  Similarly,  the  topologies  of  6x2x2  (Z- 
spacing  =  1.5  m)  and  6x2x3  (Z-spacing  =  0.75  m)  produced  similar  RMS  errors  of  3  m, 
7  m  and  11m  for  Xs,  Ys  and  Zs  readings,  respectively.  Thus,  it  is  recommended  that  two 


sensors  per  soldier  will  be  sufficient. 


□  6x2x2  (Z-spacing=1  m) 

■  6x2x3  (Z-spacing=0.5  m) 

□  6x2x2  (Z-spacing  =  1.5  m) 

□  6x2x3  (Z-spacing=  0.75  m) 


Figure  26.  RMS  error  ^for  two-sensor  (6x2x2  with  Z-spacing  of  1  m  and  1.5  m) 
and  three-sensor  (6x2x3  with  Z-spacing  of  0.5  m  and  0.75  m)  topologies, 
representing  wearable  configuration.  Results  are  obtained  by  averaging  output 
from  100  Monte  Carlo  simulations. 


4.  Effect  of  Noise 

Measurement  noise  is  expected  to  have  a  significant  effect  on  the  estimation 
accuracy  [14].  The  variance  in  acoustic  travel  time  cr,^  along  a  ray  path  from  source  to 
receiver  is  proportional  to  the  propagation  distance  as  given  by  [3] 

2  ^  2  R  L, 

=0.4z  ^ 
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where  L  represents  the  correlation  length  of  atmospheric  turbulence,  v  is  the  speed  of 
sound  and  ;^^is  the  mean-squared  sound-speed  fluctuations  near  the  ground.  However, 
the  proportionality  property  is  valid  only  for  short  distances  of  less  than  60  m.  A  higher 
order  approximation  is  required  to  provide  a  more  realistic  estimation.  In  this  study,  the 
RDOA  noise  variance  is  modeled  as  a  second-order  approximation 

where  a  is  a  constant  that  is  set  to  0.01  for  the  study. 

To  study  the  effect  of  noise  on  the  location  estimation  performance,  a  scenario 
consisting  of  24  sensor  nodes  (6x2x2),  as  shown  in  Figure  27,  with  a  10-m  spacing 
between  adjacent  sensors  is  simulated.  The  source  is  placed  100  m  from  the  origin. 


source 


(0,100,0) 


Figure  27.  Simulated  scenario  of  24  sensor  nodes  (6x2x2)  with  a  10-m  spacing 
between  adjacent  sensors  is  simulated.  The  source  is  placed  100  m  from  the 

origin. 


Figure  28  and  Table  5  show  the  error  mean  jn  ,  error  variance  (7^ ,  and  RMS  error 
Y  for  different  noise  levels. 
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method  with  different  noise  variance  <7^  for:  (a)  error  mean  // ,  (b)  error  variance 
(7^ ,  and  (c)  RMS  error  y,  representing  a  6x2x2  topology.  Results  are  based  on 

100  Monte  Carlo  runs. 
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Noise 

variance 

0.001 

0.01 

0.1 

1 

Error  mean 
/n^j  in  meters 

Xs SI 

-2.5658 

-2.5682 

-2.5822 

-1.021 

Ys SI 

11.133 

11.141 

11.115 

-2.0383 

Zs Sl 

-2.1288 

-2.1309 

-2.1485 

-1.3407 

Rs SI 

11.183 

11.191 

11.168 

-1.8955 

Bearing SI 

1.7183 

1.7198 

1.7261 

2.8112 

Error  mean 

in  meters 

Xs ML 

-2.1565 

-1.4044 

-1.4486 

0.007947 

Ys ML 

10.818 

9.0713 

9.2698 

-0.37713 

Zs ML 

-1.6793 

-0.95559 

-0.97037 

-0.01764 

Rs ML 

10.854 

9.0848 

9.2844 

-0.32973 

Bearing ML 

1.4118 

0.89127 

0.91088 

1.535 

Error  variance 
ah  in  square 
meters 

Xs SI 

2.57E-05 

0.002569 

0.25342 

11.914 

Ys SI 

0.000386 

0.038694 

3.8292 

359.54 

Zs SI 

1.84E-05 

0.001845 

0.18198 

7.5425 

Rs SI 

0.000394 

0.039445 

3.9012 

353.59 

Bearing SI 

9.72E-06 

0.000973 

0.09733 

11.599 

Error  variance 

in 

square  meters 

Xs ML 

2.66E-01 

0.04107 

0.1096 

5.4121 

Ys ML 

1.837 

1.1165 

3.4264 

196.71 

Zs ML 

2.88E-01 

0.02666 

0.063778 

2.828 

Rs ML 

1.8589 

1.1228 

3.4485 

195.4 

Bearing ML 

1.39E-01 

0.017059 

0.041712 

1.4399 

RMS  error 

Ygj  in  meters 

Xs SI 

2.5658 

2.5687 

2.6308 

3.5996 

Ys SI 

11.133 

11.143 

11.286 

19.071 

Zs SI 

2.1288 

2.1313 

2.1905 

3.0561 

Rs Sl 

11.183 

11.193 

11.341 

18.899 

Bearing SI 

1.7183 

1.7201 

1.7541 

4.4161 

RMS  error 
y^i^  in  meters 

Xs ML 

2.2173 

1.419 

1.4859 

2.3264 

Ys ML 

10.902 

9.1327 

9.4528 

14.03 

Zs ML 

1.7629 

0.96944 

1.0027 

1.6818 

Rs ML 

10.939 

9.1464 

9.4683 

13.983 

Bearing_ML 

1.4603 

0.90079 

0.9335 

1.9484 

Table  5.  Effect  of  noise  on  the  location  estimation  performance  of  the  hybrid  SI/ML 
method  with  different  noise  variance  cr^  for:  error  mean  /j. ,  error  variance  (7^ , 


and  RMS  error  y ,  representing  a  6x2x2  topology.  Results  are  obtained  by 
averaging  output  from  100  Monte  Carlo  simulations. 


It  can  be  seen  that  as  the  noise  variance  (J^  increases,  the  error  variance  (7^ 

increases.  For  example,  the  error  variance  (7^  for  the  Rs  value  at  the  ML  stage  increased 

2  2  2 

from  1.1m  to  3.4  m  when  the  noise  variance  increased  from  0.1  to  1. 
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Although  it  is  expected  that  a  noisier  environment  will  lead  to  greater  estimation 
errors,  it  is  interesting  to  note  that  at  high  noise  levels  (i.e.,  (7^  =  1 ),  ML  actually  produces 

an  error  mean  that  is  better  than  for  low  noise  levels  (i.e.,  <7^  =  0.001,  0.01  or  0.1).  For 

example,  the  error  mean  for  the  Rs  values  for  the  ML  stage  are  11m  and  -2  m  for 

values  of  0.1  and  1,  respectively.  This  anomaly  is  likely  due  to  the  inaccurate  SI  value 
that  is  used  as  the  initial  value  for  the  ML  stage.  At  low  noise  levels,  the  SI  stage 
introduces  large  bias  of  approximately  11  m  and  leads  to  similarly  large  errors  of 
approximately  10  m  for  the  ML  stage.  The  smaller  bias  of  2  m  introduced  by  SI  at  high 
noise  levels  actually  leads  to  better  estimation  accuracy  in  the  ML  stage. 

5.  Effect  of  Source-Sensor  Distance 

To  study  the  effect  of  the  source-sensor  distance  on  the  location  estimation 
performance,  a  scenario  consisting  of  24  sensor  nodes  (6x2x2)  is  simulated.  The 
distance  between  the  source  and  the  sensor  is  varied  from  10  m  to  200  m.  Figure  29  and 
Table  6  show  the  error  mean//,  error  variance  (7^,  and  RMS  error  y  for  the  various 
source-sensor  distances  in  this  range. 

As  the  source-sensor  distance  increases,  the  error  variance  (7^  increases.  It  is  seen 
that  the  error  variance  of  the  SI  method  in  the  first  stage  increases  much  rapidly  than  that 
in  the  ML  stage.  For  example,  at  the  first  SI  stage,  the  error  variance  of  Rs  value 
increases  from  4  m  to  354  m  when  the  source-sensor  distance  is  increased  from  50  m  to 
100  m;  the  error  variance  at  the  ML  stage  increases  from  2  m  to  195  m  , 
correspondingly. 

However,  when  the  source-sensor  distance  increases  from  150  m  and  beyond,  the 
ML  method  produces  very  divergent  results  due  to  the  erroneous  estimates  provided  by 
the  SI  method  and  are  not  shown  in  Figure  28. 
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Figure  29.  Effect  of  source-sensor  distance  on  location  estimation  performance  of  the 
hybrid  SI/ML  method:  (a)  error  mean  ,  (b)  error  variance  (J^ ,  and  (c)  RMS 
error  y ,  representing  a  6x2x2  topology.  Results  are  obtained  by  averaging 
output  from  100  Monte  Carlo  simulations.  Readings  for  150  m  and  beyond  are  not 

shown  due  to  erroneous  results. 
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Source¬ 

sensor 

distance 

10  m 

50  m 

100  m 

150  m 

200  m 

Error  mean 

in  meters 

Xs SI 

0.55436 

-1.9831 

-1.021 

3.6574 

7.7222 

Ys SI 

-2.3537 

3.4415 

-2.0383 

-48.569 

-112.68 

Zs Sl 

-1.6988 

-1.8335 

-1.3407 

1.7466 

4.6552 

Rs SI 

-2.1474 

3.5276 

-1.8955 

-45.186 

-99.574 

Bearing_ 

SI 

13.165 

2.8987 

2.8112 

16.813 

42.603 

Error  mean 

in  meters 

Xs ML 

-0.01412 

-0.54358 

0.007947 

1.21E-I-08 

2.76E-I-07 

Ys ML 

0.005905 

1.7556 

-0.37713 

3.36E-I-08 

4.72E-I-08 

Zs ML 

-0.0134 

-0.38831 

-0.01764 

-7.11E-I-07 

-1.56E-r08 

Rs ML 

0.005981 

1.7648 

-0.32973 

4.27E-I-08 

5.73E-I-08 

Bearing_ 

ML 

0.15755 

0.93126 

1.535 

5.4505 

10.156 

Error  variance 
in  square 

meters 

Xs Sl 

0.004201 

1.2665 

11.914 

36.192 

47.508 

Ys SI 

0.005741 

3.5715 

359.54 

3341.5 

7866.5 

Zs SI 

0.001391 

0.90292 

7.5425 

18.249 

23.135 

Rs Sl 

0.005977 

3.7969 

353.59 

2707.9 

5533.6 

Bearing_ 

SI 

0.000324 

1.9415 

11.599 

1240 

3297.6 

Error  variance 

in 

square  meters 

Xs ML 

0.000607 

0.32603 

5.4121 

2.18E-I-18 

5.71E-I-17 

Ys ML 

0.000151 

2.001 

196.71 

3.46E-I-18 

1.32E-I-18 

Zs ML 

0.000558 

0.21051 

2.828 

1.70E-I-17 

1.96E-I-17 

Rs ML 

0.000149 

2.0198 

195.4 

5.75E-I-18 

2.00E-I-18 

Bearing_ 

ML 

0.025758 

0.28883 

1.4399 

86.278 

144.19 

RMS  error 

Y^i  in  meters 

Xs SI 

0.55813 

2.2802 

3.5996 

7.0404 

10.351 

Ys SI 

2.3549 

3.9262 

19.071 

75.501 

143.4 

Zs SI 

1.6992 

2.0651 

3.0561 

4.6151 

6.6938 

Rs SI 

2.1487 

4.03 

18.899 

68.918 

124.29 

Bearing_ 

SI 

13.165 

3.2162 

4.4161 

39.021 

71.503 

RMS  error 
in  meters 

Xs ML 

0.0284 

0.78836 

2.3264 

1.48E-I-09 

7.56E-I-08 

Ys ML 

0.013616 

2.2546 

14.03 

1.89E-r09 

1.24E-I-09 

Zs ML 

0.027152 

0.60108 

1.6818 

4.18E-I-08 

4.70E-I-08 

Rs ML 

0.013599 

2.2659 

13.983 

2.44E-I-09 

1.53E-r09 

Bearing_ 

ML 

0.2249 

1.0752 

1.9484 

10.77 

15.727 

Table  6. 


Effect  o 


source-sensor  distance  on  location  estimation  performance  of  the  hybrid 


SEME  method:  error  mean  // ,  error  variance  ,  and  RMS  error  f  ^  representing 


a  6x2x2  topology.  Results  are  obtained  by  averaging  output  from  100  Monte 

Carlo  simulations. 
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Another  scenario  is  carried  out  to  simulate  the  effect  of  source-sensor  distance  on 
the  location  estimation  performance  for  a  group  of  soldiers  equipped  with  two  sensors 
each.  A  scenario  consisting  of  24  sensors  (6x2x2)  is  simulated  to  represent  twelve 
soldiers  who  are  individually  equipped  with  two  sensors  that  are  spaced  1  m  apart  on 
their  body.  Each  soldier  is  stationed  10  m  apart,  and  the  source  is  situated  50  m  to  200  m 
from  the  origin.  Figure  30  shows  the  RMS  error  y  for  tho  6x2x2  soldier  configuration 
with  different  source-sensor  distances:  50  m,  100  m,  150  m  and  200  m. 

It  is  seen  that  when  the  source-sensor  distance  increases,  the  RMS  error  y 
increases.  For  example,  at  a  source-sensor  distance  of  100  m,  the  RMS  error  yfor  Xs,  Ys 
and  Zs  readings  at  the  MF  stage  is  3  m,  8  m  and  15  m,  respectively;  at  a  source-sensor 
distance  of  200  m,  the  RMS  error  yfor  Xs,  Ys  and  Zs  readings  at  the  MF  stage  is  8  m,  45 
m  and  32  m,  respectively. 
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Figure  30.  RMS  error  ^for  different  values  of  source-sensor  distances:  50  m,  100  m, 
150  m  and  200  m,  representing  a  6x2x2  wearable  configuration.  Results  are 
obtained  by  averaging  output  from  100  Monte  Carlo  simulations. 


6.  Effect  of  Node  F allures 

Node  failures  can  occur  when  a  sensor’s  battery  runs  out  or  the  sensor  may  be 
asleep  when  the  source  is  detected.  The  scenario  representing  twelve  soldiers  with  two 
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sensors  each  (6x2x2)  is  studied.  Each  soldier  is  stationed  10  m  apart,  and  the  source  is 
100  m  from  the  origin.  Figure  31  shows  the  RMS  error  as  a  function  of  the  number  of 
sensor  node  failures. 

The  RMS  error  7  steadily  increases  as  the  number  of  sensor  node  failures 
increases  from  two  to  twelve.  For  example,  the  RMS  error  for  the  Rs  value  increases 
gradually  from  10  m  to  20  m  when  the  number  of  sensor  node  failures  increases  from 
zero  to  ten.  This  demonstrates  the  advantage  of  having  distributed  wireless  sensor 
networks  where  sensor  node  failures  lead  to  a  graceful  degradation  of  performance  and 
not  a  sudden  system  failure.  Note  that  drastic  increase  in  the  RMS  error  y  at  the  SI  stage 
is  observed  only  when  there  are  twelve  or  more  sensor  node  failures  (i.e.,  50  %  sensor 
node  failure).  At  such  a  high  failure  rate,  the  SI  method  used  in  the  first  stage  of  the 
hybrid  SFMF  method  gives  a  large  RMS  error  of  50  m  for  the  Rs  value,  but  the  MF 
method  in  the  second  stage  is  able  to  tolerate  the  errors  and  produced  RMS  errors  for  the 
Rs  value  of  less  than  10  m. 
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‘o' 


‘o' 


□  0  failure 

■  2  failures 

□  4  failures 

□  6  failures 

■  8  failures 

□  10  failures 

■  12  failures 


Figure  31.  RMS  error  y  as  the  function  of  the  number  of  sensor  node  failures, 
representing  a  6x2x2  wearable  configuration.  Results  are  obtained  by  averaging 
output  from  100  Monte  Carlo  simulations. 
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7.  Effect  of  Location  Perturbations 

In  field  deployment,  it  is  expected  that  there  will  be  some  perturbation  (or  slight 
deviation)  in  the  sensor  locations  as  perfectly  straight  line  troop  formation  or  exact 
symmetry  is  usually  not  possible.  This  perturbation  effect  on  location  estimation 
performance  is  studied  by  introducing  some  degree  of  randomness  to  the  sensor 
positions.  The  standard  deviation  of  the  random  distance  between  soldiers  is  varied  from 
0.01  to  1;  the  standard  deviation  of  the  two  sensor  locations  on  a  soldier  is  set  to  0.01. 
Figure  32  shows  the  RMS  error  ^'for  different  values  of  standard  deviations  of  location 
perturbations. 


It  is  seen  that  RMS  error  ^'is  similar  for  low  and  high  standard  deviations  of 


location  perturbation.  For  example,  the  RMS  error  for  the  Rs  value  is  between  7  m  to  9 
m  for  the  ML  method  at  the  second  stage  of  the  hybrid  SFML  algorithm  for  the  various 
levels  of  location  perturbation  standard  deviations.  This  shows  that  location  perturbation 
has  limited  effect  on  the  estimation  accuracy.  This  clearly  illustrates  the  flexibility  of  the 
network  configuration  and  demonstrates  the  suitability  for  potential  field  deployment. 


Location  Perturbation  (Standard  Deviation) 


□  stdev  =  0 

■  stdev=  0.01 

□  stdev  =  0.1 

□  stdev  =  1 


Figure  32.  RMS  error  ^'for  different  values  of  standard  deviations  of  location 
perturbations,  representing  a  6x2x2  wearable  configuration 
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C.  WIRELESS  NETWORK  SIMULATION 

The  wireless  network  simulation  is  carried  out  using  NS-2  simulator,  a  discrete 
event  simulator  suitable  for  networking  research.  The  wireless  extension  and  the  802.15.4 
extension  [19]  are  used  to  conduct  the  studies  on  scenarios  using  an  IEEE  802.1 1-based 
wireless  local  area  network  (WEAN)  and  an  IEEE  802.15.4-based  wireless  personal  area 
network  (WPAN).  Eigure  33  shows  the  protocol  architecture  of  the  NS-2  wireless 
extension  to  support  WPAN  simulation. 


Wireless  Scenario 
Definition 


NS-2 


Eigure  33.  The  protocol  architecture  of  the  IEEE  802. 15.4  WPAN  extension  in  NS-2 

(from  [8]). 


A  scenario  of  a  wireless  sensor  network  consisting  of  26  sensor  nodes  situated 
over  a  50  mx50  m  area  is  simulated  using  NS-2.  Eigure  34  shows  the  layout  of  the 
sensors  defined  for  the  scenario  used  for  the  simulation.  The  traffic  scenario  of  twelve 
sensor  nodes  (#1  to  #12)  sending  Constant  Bit  Rate  (CBR)  traffic  to  a  central  node  (#19) 
at  the  same  time  is  simulated.  The  Adhoc  On-demand  Distance  Vector  (AODV)  routing 
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protocol  [20]  is  run  on  top  of  the  802.15.4  MAC  and  PHY  layers  to  study  the 
performanee  of  WPAN  networks.  Similarly,  the  AODV  routing  protoeol  is  run  on  top  of 


802. 1 1  MAC  and  PHY  layer  to  study  the  performance  of  WLAN  networks. 


nam:  wpcin_cbr.nam 

1  RIe  ^ews  Analysis 

wpan_cbr.nain 

1 

■ 

!  >■  1 

60  002000 

Step;  100  Oms 

1 1 .0)  cbr  traffic  from  node  #1  to  12  to  node  #19 


Figure  34.  NS-2  WPAN  and  WLAN  Simulation  Seenario:  26  sensor  nodes  in  a  50 
x50  m  area.  The  Time  Interval  Length  (TIL)  used  in  the  simulation  is  equal  to 

one  second. 


1.  Sending/Receiving  Throughput 

Figure  35  shows  the  sending  and  receiving  throughput  for  WPAN  and  WLAN 
networks.  The  Time  Interval  Length  (TIL)  used  in  the  simulation  is  equal  to  one  second. 
The  average  sending  throughput  refers  to  the  average  throughput  of  nodes  1  to  12  that  are 
sending  CBR  traffic  to  central  node  19;  the  average  receiving  throughput  refers  to  the 
average  throughout  received  by  node  19. 

For  WPAN  networks,  the  average  sending  throughput  is  approximately  100  kbps 
and  the  receiving  throughput  is  approximately  20  kbps.  It  is  observed  that  the  WLAN 
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network  has  better  throughput  than  WPAN  network.  For  WLAN  networks,  the  average 
sending  throughput  is  approximately  300  kbps  and  the  receiving  throughput  is 
approximately  200  kbps. 


x10° 


(a) 


X  10*' 


(b) 

Figure  35.  Sending  and  Receiving  Throughput  of:  (a)  WPAN  and  (b)  WLAN.  The 
Time  Interval  Length  (TIL)  used  in  the  simulation  is  equal  to  one  second.  NS-2 
simulation  of  a  26-node  network  operating  over  an  area  of  50  x  50  m^. 
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2.  End-to-End  Delay 

Figure  36  shows  the  histograms  and  the  summary  statistics  for  the  end-to-end 
delay  for  the  WPAN  and  WLAN  networks.  For  WPAN  networks,  approximately  90  %  of 
the  traffic  has  end-to-end  delay  of  less  than  10  ms.  The  average  WPAN  end-to-end  delay 
is  6.4  ms.  The  minimum  WPAN  delay  is  0.35  ms  and  the  maximum  delay  is  1 15  ms. 

It  is  observed  that  the  end-to-end  delay  in  the  WLAN  network  is  slightly  higher 
than  that  in  the  WPAN  network.  For  WLAN  networks,  approximately  97  %  of  the  traffic 
has  end-to-end  delay  of  less  than  50  ms.  The  average  WLAN  end-to-end  delay  is  17  ms. 
The  minimum  WLAN  delay  is  0.2  ms  and  the  maximum  delay  is  564  ms.  This 
observation  of  lower  end-to-end  delay  in  the  WPAN  networks  compared  to  the  WLAN 
networks  is  consistent  with  the  results  of  the  WPAN  performance  study  performed  by 
Zheng  [19]. 


Simulation  End2End  delays  in  seconds: 


Minimal  delay  (CN,0N,PID):  0.000352009  (5,1 9.928] 
Maximal  delay  (CN.ON.PID):  0.115648059  (7,19,2850) 
Average  delay:  0.006437120238 


(a) 
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Simulation  End2End  delays  in  seconds: 


Minimal  delay  (CN.ONPID):  0.000230033  (6,-1 ,0) 
Maximal  delay  (CN,ON,PID):  0.564433939  (10,19,93) 
Average  delay;  0.01735257521 


(b) 

Figure  36.  End-to-End  delay  for:  (a)  WPAN  and  (b)  WEAN.  NS-2  simulation  of  a 
26-node  network  operating  over  an  area  of  50  x  50  m  . 


3.  Packet  Size 

Eigure  37  shows  the  simulation  summary  statistics  for  WPAN  and  WEAN 
networks.  The  average  WPAN  packet  size  is  51.37  bytes  with  a  minimum  packet  size  of 
5  bytes  and  a  maximum  packet  size  of  111  bytes.  The  WEAN  networks  have  a  similar 
packet  size  as  WPAN  networks.  The  average  WEAN  packet  size  is  57.51  bytes  with  a 
minimum  packet  size  of  28  bytes  and  a  maximum  packet  size  of  142  bytes. 
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Simulation  information; 

Simulation  information; 

Simulation  length  in  seconds: 

58.99960653 

Simulation  length  in  seconds; 

58.86696361 

Number  of  nodes: 

26 

Number  of  nodes: 

26 

Number  of  sending  nodes: 

26 

Number  of  sending  nodes: 

26 

Number  of  receiving  nodes; 

23 

Number  of  receiving  nodes: 

23 

Number  of  generated  packets; 

23206 

Number  of  generated  packets: 

42689 

Number  of  sent  packets: 

21326 

Number  of  sent  packets; 

42250 

Number  of  forwarded  packets: 

0 

Number  of  forwarded  packets: 

0 

Number  of  dropped  packets; 

16527 

Number  of  dropped  packets; 

7929 

Number  of  lost  packets: 

6866 

Number  of  lost  packets; 

11339 

Minimal  packet  size; 

5 

Minimal  packet  size; 

28 

Maximal  packet  size; 

111 

Maximal  packet  size; 

142 

Average  packet  size; 

51.3744 

Average  packet  size: 

57.5045 

Number  of  sent  bytes: 

1151295 

Number  of  sent  bytes: 

2635686 

Number  of  forwarded  bytes; 

0 

Number  of  forwarded  bytes: 

0 

Number  of  dropped  bytes; 

1059966 

Number  of  dropped  bytes: 

478598 

Packets  dropping  nodes: 

|01  23456789  10 

Packets  dropping  nodes: 

|01  2345789  10  r 

(a)  (b) 

Figure  37.  Simulation  summary  statistics  for:  (a)  WPAN  and  (b)  WLAN  networks. 
NS-2  simulation  of  a  26-node  network  operating  over  an  area  of  50  x  50  m^. 

4.  Jitter 

Figure  38  shows  the  jitter  plots  of  transmitted  and  received  traffic  for  WPAN  and 
WLAN  networks.  It  is  observed  that  the  packet  jitter  is  low  for  both  WLAN  and  WPAN 
networks. 

For  the  WPAN,  the  maximum  jitter  of  sending  packets  is  200  ms  and  the 
maximum  jitter  of  receiving  packets  is  320  ms.  The  average  jitter  of  sending  and 
receiving  packets  is  approximately  8  ms.  For  WLAN,  the  maximum  jitter  of  sending 
packets  is  200  ms  and  the  maximum  jitter  of  receiving  packets  is  also  200  ms.  The 
average  jitter  of  sending  and  receiving  packets  is  approximately  5  ms. 


61 


jitter  of  received  packet  [sec] 


(a) 
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0  50  100  150  200  250 


sequence  number 


(c) 


0  50  1M  150  200  250 


sequence  number 


(d) 

Figure  38.  Jitter  plot  of:  (a)  all  sent  WPAN  packets,  (b)  all  received  WPAN  packets 
(c)  all  sent  WLAN  packets,  (d)  all  received  WLAN  packets.  NS-2  simulation  of  a 
26-node  network  operating  over  an  area  of  50  x  50  m^. 


The  simulation  results  show  that  the  wireless  sensor  network  has  a  throughput 

between  20  kbps  and  100  kbps,  a  low  average  delay  of  less  than  10  ms  and  a  low  average 

jitter  in  the  range  of  5  ms  to  10  ms.  Assuming  that  the  payload  sent  by  each  sensor  to  a 
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server  for  proeessing  is  1  kByte,  as  has  been  implemented  by  Duekworth  et  al.  [6],  the 
required  transmission  time  is  approximately  between  0.08  s  and  0.4  s.  After  considering 
the  processing  time  and  other  delays,  the  total  delay  can  be  approximated  to  be  less  than 
one  second.  This  indicates  that  the  wireless  sensor  network  is  suitable  for  implementation 
in  fast-response  counter-sniper  applications. 

This  chapter  has  presented  the  simulation  results  of  the  hybrid  SI/ML  estimation 
method  used  in  the  two-step  localization  process  for  counter- sniper  applications. 
Simulation  results  indicate  that  the  proposed  two-step  localization  method  is  suitable  for 
counter-sniper  applications.  For  example,  twelve  soldiers  with  two  sensors  each  (6x2x2) 
can  detect  a  sniper  100  m  away  with  an  RMS  error  of  3  m,  8  m  and  15  m  for  the  Xs,  Ys 
and  Zs  values,  respectively.  Better  accuracy  can  be  obtained  by  increasing  the  number  of 
sensors  or  by  increasing  the  inter-sensor  spacing.  The  NS-2  based  simulation  results 
indicate  that  both  WPAN  and  WLAN  networks  are  suitable  for  implementation  in  fast- 
response  counter-sniper  applications  due  to  their  high  throughput,  low  delay  and  low 
jitter.  The  next  chapter  will  conclude  the  thesis  by  summarizing  the  significant  results  and 
proposing  possible  future  work. 
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V.  CONCLUSION 


This  thesis  focused  on  source  localization  using  wireless  sensor  networks. 
Wireless  sensor  nodes  worn  on  soldiers  or  installed  on  vehicles  can  be  deployed  to  form  a 
distributed  sensor  array  to  locate  the  source  of  sniper  fire.  A  two-step  source  localization 
process  is  proposed  for  counter-sniper  applications:  the  time  difference  of  arrival 
(TDOA)  values  are  first  determined  using  the  generalized  cross  correlation  (GCC) 
method  and  then  the  TDOA  values  are  used  by  a  hybrid  spherical  interpolation/maximum 
likelihood  (SI/ML)  estimation  method  to  determine  the  shooter  location.  A  simulation 
model  has  been  developed  in  MATLAB  to  study  the  performance  of  the  hybrid  SI/ML 
estimation  method.  A  wireless  sensor  network  is  simulated  in  NS-2  to  study  the  network 
throughput,  delay  and  jitter  of  the  IEEE  WPAN  and  WEAN  networks. 

A.  SIGNIFICANT  RESULTS 

The  effects  of  number  of  sensor  nodes,  inter-sensor  spacing,  topology,  noise, 
source-sensor  distance,  sensor  node  failure,  and  location  perturbation  on  the  location 
estimation  performance  of  the  hybrid  SI/ME  method  have  been  studied  through 
MATEAB  simulation. 

Simulation  results  show  that  the  error  variance  can  be  reduced  by  increasing  the 
number  of  sensor  nodes  or  the  inter-sensor  spacing.  To  increase  the  estimation  accuracy, 
the  number  of  sensors  as  well  as  the  inter-sensor  spacing  should  be  increased. 

Eor  the  case  of  wearable  sensors,  the  constraint  of  small  inter-sensor  spacing  on 
the  human  body  degrades  the  estimation  accuracy.  On  the  other  hand,  vehicular 
configuration  providing  larger  vertical  and  horizontal  inter-sensor  spacing  can  help 
improve  the  estimation  accuracy.  We  see  that  the  configuration  of  a  group  of  twelve 
soldiers  or  vehicles,  with  two  wireless  sensors  per  soldier  or  vehicle,  is  sufficient  to 
accurately  detect  a  sniper  situated  100  m  away  with  an  RMS  error  of  3  m,  8  m  and  15  m 
for  the  Xs,  Ys  and  Zs  values,  respectively.  By  using  a  longer  range  wireless  network  such 
as  WEAN,  the  inter-sensor  spacing  can  be  increased  and  provides  better  estimate 
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accuracy  to  detect  a  sniper  located  farther.  Simulation  results  also  indicate  that  the  sensor 
topology  should  be  well  represented  in  all  three  dimensions  to  obtain  desirable  resolution. 

The  source  estimation  accuracy  is  not  adversely  affected  when  encountering  node 
failures  or  location  perturbations.  Results  show  that  sensor  node  failures  lead  to  a  gradual 
degradation  of  accuracy  and  not  a  sudden  total  system  failure.  This  demonstrates  the 
advantage  of  using  distributed  wireless  sensor  networks  to  provide  greater  system 
robustness.  It  is  observed  that  location  perturbations  do  not  adversely  affect  the  location 
estimation.  This  illustrates  the  flexibility  of  the  network  configuration  and  demonstrates 
the  suitability  for  potential  field  deployment. 

Network  simulation  results  based  on  the  NS-2  package  show  that  the  wireless 
sensor  network  can  support  a  throughput  of  between  20  kbps  and  100  kbps,  a  low  average 
delay  of  less  than  10  ms,  and  a  low  average  jitter  in  the  range  of  5  ms  to  10  ms.  These 
are  suitable  for  fast  response  network-based  sniper  detection  applications  that  require  low 
delay  and  fast  data  exchange. 

B.  FUTURE  WORK 

This  thesis  dealt  with  a  single  source  scenario.  In  a  future  effort,  the  work 
reported  here  can  be  extended  to  multiple  sources.  Since  the  SI  method  is  only  applicable 
to  a  single  source  case,  a  different  estimation  method  should  be  considered.  Other  source 
types,  such  as  RF  and  seismic  signals,  can  also  be  included  in  the  simulation  model. 

We  have  seen  that  the  SI  method  provided  good  initial  values  for  the  ML  method. 
However,  the  SI  method  can  introduce  bias  at  low  noise  levels  and  affect  the  estimation 
accuracy  of  the  ML  method.  Alternative  estimation  methods,  such  as  linear  interpolation 
(LI)  [21]  and  cluster-based  schemes  [16],  can  be  studied  to  replace  the  SI  method  in  the 
two-step  process.  We  also  observed  that  the  ML  method  failed  to  converge  when  the  SI 
method  produced  noisy  location  estimates.  Further  study  is  required  to  analyze  this 
behavior. 

To  improve  the  estimation  accuracy  of  a  nearby  source,  beamforming  techniques 
can  be  applied  as  they  generally  provide  a  more  accurate  and  robust  performance  in  the 
near-field  [17] [22]  [23].  As  discussed  in  this  thesis,  the  proposed  two-step  process  is 
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suitable  for  far- field  applieations.  A  eombination  of  beamforming  and  the  two-step 
process  can  be  studied  in  a  future  effort  to  cover  both  near- field  and  far- field  sources. 

In  this  work,  the  two-step  source  localization  process  was  simulated  in  MATLAB. 
This  work  can  be  extended  by  implementing  the  scheme  in  hardware  using  Crossbow 
motes  and  conducting  field  experiments.  The  measured  results  can  be  used  to  validate  the 
simulated  results  provided  in  this  thesis. 
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APPENDIX:  MATLAB  SOURCE  CODE 


This  appendix  lists  the  MATLAB  source  code  used  to  study  the  performance  of 
the  Spherical  Interpolation  (SI)  and  Maximum  Likelihood  (ML)  source  localization 
estimation  method. 

% 

%  Calculate  Source  Location  based  on  Range  Difference  of  Arrival  (RDOA) 

%  (1)  Least  Square  Estimation  -  Spherical  Interpolation  (SI)  Method 
%  (2)  Maximum  Likelihood  (ML)  Method 

% 

%  N  sensors,  1  source 

%  Using  sensor  1  as  reference  i.e  (xl=0,  yl=0,  zl=0) 

% 

% 

clear; 

% 

%  Definition 

% 

nRun=100;  %  number  of  Monte  Carlo  runs 

%  uncomment  one  of  them 
%  bML=0;  %  turn  off  ML  calculation 

bML=l;  %  turn  on  ML  calculation 

%  uncomment  one  of  them 

perturb=0;  %  turn  off  location  perturbation 

%  perturb=l;  %  turn  on  location  perturbation 


% 

%  Actual  source  location  (m)  in  Cartesian  coordinates  x,  y  and  z 
%  Note:  For  simplicity,  we  only  varies  y  for  our  simulation 

% 

xs_src_actual=[0] ; 

%  Varies  the  Y  position  (Choose  1) 

% 

ys_src_actual=[100];  %100  m 
zs_src_actual=[0] ; 

Rs_actual=sqrt(xs_src_actual.^2  +  ys_src_actual.^2  +  zs_src_actual.^2); 
%  calculate  corresponding  range  Rs 

bearing_actual=[xs_src_actual;  ys_src_actual;  xs_src_actual]/Rs_actual; 
%  calculate  corresponding  bearing 
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% 

%  Actual  sensor  location  (m)  in  Cartesian  coordinates  x,  y  and  z 
%  Note:  For  simplicity,  we  only  use  integers  and  then  multiply  with 
%  a  scaling  factor  to  produce  the  actual  coordinates. 

%e.g.  [5  10  15]  =  [1  2  3]  *5  ; 

% 


%  Scale  wrt  to  Im  (Choose  1) 

% - 

scale_dist  =  10  %  10  m 

%  (Choose  1  of  the  following  sensor  configuration  for  study) 

% 

%  12x2  sensors  arranged  2  rows 

%  xi=[0  01122334455  001122334455  ].*scale_dist; 

%  yi=[0  10101010101  01010101010  l].*scale_dist; 

%  zi=[0  0  0  0  0  0  0  0  0  0  0  0  1  1  1  1  1  1  1  1  1  1  1  1].*1.0; 

%  Soldier  configuration  (each  with  2  sensors).  z=  Im  apart  vertically 

temp=size(xi); 

nSen=temp(l,2);  %  number  of  sensor  (>4) 

noisestd=l; 
if  (perturb==l) 
randn('state',0); 
tmpl=randn(3,  nSen); 
for  i=l:2:nSen 

xi(i)=  xi(i)  +  noisestd*tmpl(l,i); 

xi(i+l)=  xi(i+l)+(noisestd+0.01)*tmpl(l,i);  %  less  variance  on  the  body 
yi(i)=  yi(i)  +  noisestd*tmpl(2,i); 

yi(i+l)=  yi(i+l)+(noisestd+0.01)*tmpl(2,i);  %  less  variance  on  the  body 
end 

zi=zi+0.01*tmpl(3,:);  %  less  variance  on  the  body 
end 

%  RD  noise  (Choose  1) 

% 

Noise_Factor=0.001;  %  noise  std  =  Std_Norm  *  (source  distance).  % 
we  expect  bigger  noise  variance  for  larger  distance. 
Noise_Var=(Noise_Factor*Rs_actual)^2; 

% 

%  Functions 

% 
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%  Random  Process 
%  AWGN 
randn('state',0); 

noise  =  sqrt(Noise_Var)*randn(nRun,  1); 

%noise_mean  =  mean(noise,  2);  %  average  along  row 

for  k=l:nRun,  %  Monte  Carlo  Simulation 

Xi=[xi'  yi'  zi']; 

Di=  sqrt  ((xi-xs_src_actual).^2  +  (yi-ys_src_actual).^2  +  (zi-zs_src_actual).^2); 
Ri=  sqrt  ((xi).^2  +  (yi).^2  +  (zi).^2); 


locSen=[xi'  yi'  zi']; 


%  using  N  sensors 
for  i=l:nSen-l 

%d21=Di(2)-Di(l); 

%d31=Di(3)-Di(l);.. 

%dnl=Di(n)-Di(l); 

%d=[d21;d31;...;dnl]; 

d(i,l)=Di(i+l)-Di(l)+noise(k);  %add  noise  to  RD  estimates 

%  delta2=Ri(2)^2-d(l)^2; 

%  delta3=Ri(3)^2-d(l)^2;... 

%  deltan=Ri(n)^2-d(l)^2; 

%  delta=  [delta2  ;delta3 ; . . .  deltan] ; 
delta(i,l)=Ri(i+l)^2-d(l)^2; 

%  s2=  [xi(2)  yi(2)  zi(2)]; 

%  s3=  [xi(3)  yi(3)  zi(3)];... 

%  sn=  [xi(n)  yi(n)  zi(n)] ; 

%  s=[s2;s3;...sn]; 
s(i,:)=[xi(i+l)  yi(i+l)  zi(i+l)]; 

end 

%  define  weight  (positive  definite,  i.e.,  diagonal  positive  and  symmetrical) 
w=eye(nSen-l);  %  set  to  identity  matrix  for  unweighted  case 

Sw=(s'*w*s)^(-  l)*s'*w; 

Ps=s*Sw; 

Ps_ortho=eye(nSen- 1  )-Ps ; 

% 

%SI  method 

% 
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Rs_SI_cal=0.5*(d'*Ps_ortho*w*Ps_ortho*delta)/(d'*Ps_ortho*w*Ps_ortho*d); 

%  Calculate  Xs  for  SI  method 

Xs_row_SI  =  0.5*Sw*(delta-2*Rs_SI_cal*d); 

Xs_SI(k,:)=Xs_row_Sr  ; 

Rs_SI(k,:)=sqrt(Xs_SI(k,l)^2  +  Xs_SI(k,2)^2  +  Xs_SI(k,3)^2); 


bearing_SI(k,:)=Xs_SI(k,:)/Rs_SI(k,:); 

%  error_row=  delta  -  2*Rs_SI*d  -  2*s*Xs_row_SI;  %error 
%  error(k,:)=error_row'; 

%  % - 

%  %  Maximum  Likelihood  Method 
%  %  Objective  function  is  contained  in  mlobjfun.m 

%  % - 

if  (bML==l) 

xO  =  Xs_SI(k,:);  %  As  value  obtained  from  SI  as  starting  guess 

%  xO  =  [0  ys_src_actual  0];  %  Starting  guess 

options  =  optimset('LargeScaleVoff); 

%  %  LevenbergMarquardt 

options=optimset(options,'LevenbergMarquardtVon');  %  LM 
%  options=optimset(options,  'LevenbergMarquardtVoff);  %  Gauss  Newton 
[x,resnorm, residual, exitflag,output]= 

lsqnonlin(  @  nilobjfun,xO,  [] ,  [] , options  ,locSen,Noise_V  ar,d) ; 
Xs_ML(k,:)=x; 

Rs_ML(k,:)=sqrt(Xs_ML(k,l)^2+Xs_ML(k,2)^2+Xs_ML(k,3)^2); 

bearing_ML(k,:)=Xs_ML(k,:)/Rs_ML(k,:); 

end 


% 

%  Calculate  bias  (i.e.,  errors)  for  source  location,  range  and  bearing 
%  estimates 
% 

%  Calculate  mean 

% 

%  SI 

bias_Xs_SI(k,  I  )=Xs_SI(k,  I  )-xs_src_aetual; 
bias_Xs_SI(k,2)=Xs_SI(k,2)-ys_src_actual; 
bias_Xs_SI(k,3)=Xs_SI(k,3)-zs_src_actual; 


%  ML 

if  (bML==I) 
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bias_Xs_ML(k,  1  )=Xs_ML(k,  1  )-xs_src_actual; 
bias_Xs_ML(k,2)=Xs_ML(k,2)-ys_src_actual; 
bias_Xs_ML(k,  3  )=Xs_ML(k,3 )  -zs_src_actual ; 
end 
end 

bias_Rs_SI  =  Rs_SI-Rs_actual; 

bias_bearing_SI  =  180/pi*aeos(bearing_SI*bearing_aetual); 


if  (bML==l) 

bias_Rs_ML=Rs_ML-Rs_actual; 

bias_bearing_ML  =  180/pi*aeos(bearing_ML*bearing_aetual); 
end 

meanxs_SI=mean(bias_Xs_SI( :,!)); 
meany  s_SI=mean(bias_Xs_SI( :  ,2)) ; 
meanzs_SI=mean(bias_Xs_SI(:,3)); 
meanrs_SI=mean(bias_Rs_SI) ; 
meanbear_SI=mean(bias_bearing_SI); 

veet_mean_SI=[meanxs_SI;meanys_SI;meanzs_SI;meanrs_SI;meanbear_SI]; 

%ML 

if  (bML==l) 

meanxs_ML=mean(bias_Xs_ML(: ,  1 )) ; 
meanys_ML=mean(bias_Xs_ML(:  ,2)) ; 
meanzs_ML=mean(bias_Xs_ML( :  ,3)) ; 
meanrs_ML=mean(bias_Rs_ML) ; 
meanbear_ML=mean(bias_bearing_ML); 

vect_mean_ML=[meanxs_ML;meanys_ML;meanzs_ML;meanrs_ML; 

meanbear_  ML]; 
end 


%  Caleulate  Varianee  =  E[(a  -  mean)^2] 

% 

varxs_SI=var(bias_Xs_SI( :,!)); 
vary  s_SI=var(bias_Xs_SI( :  ,2)) ; 
varzs_SI=var(bias_Xs_SI( :  ,3)) ; 
varrs_SI=var(bias_Rs_SI); 
varbear_SI=var(bias_bearing_SI); 

veet_var_SI=[varxs_SI;varys_SI;varzs_SI;varrs_SI;varbear_SI]; 


%ML 

if  (bML==l) 
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varxs_ML=var(bias_Xs_ML( : ,  1 )) ; 
vary  s_ML=var(bias_Xs_ML( :  ,2)) ; 
varzs_ML=var(bias_Xs_ML(:,3)); 
varrs_ML=var(bias_Rs_ML) ; 
varbear_ML=var(bias_bearing_ML); 

vect_var_ML=[varxs_ML;varys_ML;varzs_ML;varrs_ML;varbear_ML]; 

end 

%  Calculate  second  moment  (RMS)=  sqrt  {E[a^2]}  =  sqrt  {mean^2  +  variance} 

% 

rmsxs_SI=sqrt(mean(bias_Xs_SI( : ,  1 )) .  ^2+varxs_SI) ; 

rmsys_SI=sqrt(mean(bias_Xs_SI(:,2)).^2+varys_SI); 

rmszs_SI=sqrt(mean(bias_Xs_SI(:,3))-^2+varzs_SI); 

rmsrs_SI=sqrt(mean(bias_Rs_SI).'^2+varrs_SI); 

rmsbear_SI=sqrt(mean(bias_bearing_SI).^2+varbear_SI); 

vect_rms_SI=[rmsxs_SI;rmsys_SI;rmszs_SI;rmsrs_SI;rmsbear_SI]; 

%ML 

if  (bML==l) 

rmsxs_ML=sqrt(mean(bias_Xs_ML( : ,  1 ) ) .  ^2+ varxs_ML) ; 
rmsy  s_ML=sqrt(mean(bias_Xs_ML( :  ,2)) .  ^2+vary  s_ML) ; 
rmszs_ML=sqrt(mean(bias_Xs_ML( :  ,3))  •  ^2+ varzs_ML) ; 
rmsrs_ML=sqrt(mean(bias_Rs_ML).^2+varrs_ML); 
rmsbear_ML=sqrt(mean(bias_bearing_ML).^2+varbear_ML); 

vect_rms_ML=[rmsxs_ML;rmsys_ML;rmszs_ML;rmsrs_ML;rmsbear_ML]; 

end 

%  Calculate  Cramer  Rao  Bound 

% 


cov_mat=Noise_Var.*(0.5*ones(length(d))+0.5*eye(length(d))); 
for  i=l:length(d) 

a  1 =[xs_src_actual-locSen(i+ 1,1)  y s_src_actual-locSen(i+ 1 ,2) 
zs_src_actual-locSen(i+ 1,3)]; 

a2=sqrt((xs_src_actual-locSen(i+ 1 , 1  ))^2+(ys_src_actual-locSen(i+ 1 ,2))^2 
+(zs_src_actual-locSen(i+l,3))^2); 
b  1 =[xs_src_actual-locSen(  1,1)  y s_src_actual-locSen(  1,2) 
zs_src_actual-locSen(  1,3)]; 

b2=sqrt((xs_src_actual-locSen(  1 , 1  ))^2+(ys_src_actual-locSen(  1 ,2))^2 
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+(zs_src_actual-locSen(  1 ,3))^2) ; 
jacobian(i,:)=  (al/a2)-(bl/b2); 
end 

fisher=jacobian'*inv(cov_mat)*jacobian; 
crlb=  trace(fisher^-l);  %  compare  with  MSE  of  Rs 

% 

%  Generate  Plots 

% 

%  hfigl=figure; 
if  (bML==l) 

plot(xi,  yi,'kv',  xs_src_actual,  ys_src_actual,  Xs_SI(:,l),  Xs_SI(:,2), 
'ko',Xs_ML(:,l),  Xs_ML(:,2),  'kd');  %  plot  both  SI  and  ML 
else 

plot(xi,  yi,'kv',  xs_src_actual,  ys_src_actual,  'k'^',  Xs_SI(:,l),  Xs_SI(:,2),  'ko'); 

%  plot  just  SI  only 
end 

title('Sensor  and  Source  Location'); 

strl=sprintf('[Xs,  Ys,  Zs,  Rs,  Bearing],  Noise  Std  =  %s*Rs',Noise_Lactor); 
str2=sprintf('SI  Method'); 

str3=sprintf('RMS  =  [%s,  %s,  %s,  %s,  %s]',  rmsxs_SI,  rmsys_SI,  rmszs_SI, 
rmsrs_SI,  rmsbear_SI); 

str4=sprintf('Mean  =  [%s,  %s,  %s,  %s,  %s]',  meanxs_SI,  meanys_SI,  meanzs_SI, 
meanrs_SI,  meanbear_SI); 

str5=sprintf('Variance  =  [%s,  %s,  %s,  %s,  %s]',  varxs_SI,  varys_SI,  varzs_SI, 
varrs_SI,  varbear_SI); 
if  (bML==l) 

str6=sprintf('ML  Method'); 

str7=sprintf('RMS  =  [%s,  %s,  %s,  %s,  %s]',  rmsxs_ML,  rmsys_ML,  rmszs_ML, 
rmsrs_ML,  rmsbear_ML); 

str8=sprintf('Mean  =  [%s,  %s,  %s,  %s,  %s]',  meanxs_ML,  meanys_ML, 
meanzs_ML,  meanrs_ML,  meanbear_ML); 

str9=sprintf(' Variance  =  [%s,  %s,  %s,  %s,  %s]',  varxs_ML,  varys_ML, 
varzs_ML,  varrs_ML,  varbear_ML); 

str=sprintf('%s  \n%s  \n%s  \n%s  \n%s  \n%s  \n%s  \n%s  \n%s',  strl,  str2,  str3,  str4, 
str5,  str6,  str7,  str8,  str9); 

legend('sensor  location',  'actual  source  location  ',  'calculated  source  location 
(SI)','calculated  source  location  (ML)'); 
else 

str=sprintf('%s  \n%s  \n%s  \n%s  \n%s  \n%s  \n%s  \n%s  \n%s',  strl,  str2,  str3,  str4, 
str5); 

legend('sensor  location',  'actual  source  location',  'calculated  source  location  (SI)'); 
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end 


text(scale_dist/10,ys_src_actual,str); 
xlabel('Distance  (metres)  in  X  direetion'); 
ylabel('Distanee  (metres)  in  Y  direetion'); 

%  generate  results  output  files 
fid  =  fopen('results.txt','w'); 

for  k=l:nRun, 

fprintf(fid,'%e\t%e\t%e\t%e\t%e\n',bias_Xs_SI(k,l),bias_Xs_SI(k,2), 
bias_Xs_SI(k,3),  bias_Rs_SI(k),  bias_bearing_SI(k)); 
end 

fprintf(fid,'\n%e\t  %e\t  %e\t  %e\t  %e\n',  meanxs_SI,  meanys_SI,  meanzs_SI, 
meanrs_SI,  meanbear_SI); 

fprintf(fid,'%e\t  %e\t  %e\t  %e\t  %e\n',  varxs_SI,  varys_SI,  varzs_SI,  varrs_SI, 
varbear_SI); 

fprintf(fid,'%e\t  %e\t  %e\t  %e\t  %e\n',  rmsxs_SI,  rmsys_SI,  rmszs_SI,  rmsrs_SI, 
rmsbear_SI); 

felose(fid); 
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